ANOVA với thiết kế khối và các biện pháp lặp lại

Question

Tôi đang cố gắng chạy một số phân tích thống kê về thử nghiệm thực địa được xây dựng trên 2 trang web trong cùng một mùa phát triển.

Tại cả hai địa điểm (Site, cấp: HF | NW) thiết kế thử nghiệm là một RCBD với 4 (n = 4) khối (Block, mức: 1 | 2 | 3 | 4 trong mỗi Site). Có 4 phương pháp điều trị - 3 dạng phân bón nitơ khác nhau và kiểm soát (không phân đạm) (Treatment, mức độ: AN, U, IU, C). Trong thời gian thử nghiệm thực địa, có 3 giai đoạn riêng biệt bắt đầu bằng việc bổ sung phân bón và kết thúc bằng việc thu hoạch cỏ. Những giai đoạn này đã được đưa ra mức 1 | 2 | 3 theo hệ số N_app.

Có một loạt các phép đo mà tôi muốn thử nghiệm các giả thuyết H0 sau trên:

Treatment (H0) không có hiệu lực về đo lường

Hai trong số các phép đo tôi đặc biệt quan tâm là : năng suất cỏ và khí thải amoniac.

Bắt đầu với năng suất cỏ (Dry_tonnes_ha) như hiện here, a nice balanced data set

Các dữ liệu có thể được tải về trong R sử dụng đoạn mã sau:

library(tidyverse) 

download.file('https://www.dropbox.com/s/w5ramntwdgpn0e3/HF_NW_grass_yield_data.csv?raw=1', destfile = "HF_NW_grass_yield_data.csv", method = "auto") 
raw_data <- read.csv("HF_NW_grass_yield_data.csv", stringsAsFactors = FALSE) 

HF_NW_grass <- raw_data %>% mutate_at(vars(Site, N_app, Block, Plot, Treatment), as.factor) %>% 
    mutate(Date = as.Date(Date, format = "%d/%m/%Y"), 
     Treatment = factor(Treatment, levels = c("AN", "U", "IU", "C"))) 

Tôi đã có một đi vào chạy một ANOVA về việc sử dụng này cách tiếp cận sau:

model_1 <- aov(formula = Dry_tonnes_ha ~ Treatment * N_app + Site/Block, data = HF_NW_grass, projections = TRUE) 

Tôi có một vài lo ngại về điều này.

Thứ nhất, cách tốt nhất để kiểm tra các giả định là gì? Đối với ANOVA một chiều đơn giản, tôi sẽ sử dụng shapiro.test() và bartlett.test() trên biến phụ thuộc (Dry_tonnes_ha) để đánh giá tính bình thường và tính không đồng nhất của phương sai. Tôi có thể sử dụng cùng một cách tiếp cận ở đây không?

Thứ hai, tôi lo ngại rằng N_app là phép đo lặp lại vì cùng một phép đo được lấy từ cùng một lô trên 3 giai đoạn khác nhau - cách tốt nhất để xây dựng các biện pháp lặp lại này vào mô hình là gì?

Thứ ba, tôi không chắc chắn cách tốt nhất để làm tổ Block trong phạm vi Site. Tại cả hai trang web, các cấp độ của Block là 1: 4. Tôi có cần phải có các cấp độ Block duy nhất cho mỗi trang web không?

Tôi có another data set for NH3 emissions here. Mã R để tải xuống:

download.file('https://www.dropbox.com/s/0ax16x95m2z3fb5/HF_NW_NH3_emissions.csv?raw=1', destfile = "HF_NW_NH3_emissions.csv", method = "auto") 
raw_data_1 <- read.csv("HF_NW_NH3_emissions.csv", stringsAsFactors = FALSE) 

HF_NW_NH3 <- raw_data_1 %>% mutate_at(vars(Site, N_app, Block, Plot, Treatment), as.factor) %>% 
    mutate(Treatment = factor(Treatment, levels = c("AN", "U", "IU", "C"))) 

Vì điều này, tôi có tất cả những lo ngại ở trên với việc bổ sung tập dữ liệu không cân bằng. Tại HF cho N_app 1 n = 3, nhưng đối với 3 n = 4 Tại NW n = 4 cho tất cả các mức N_app. Tại NF phép đo được chỉ thực hiện trên Treatment mức U và IU Tại NW measuremnts được thực hiện trên Treatment mức AN, U và IU

Tôi không chắc chắn làm thế nào để đối phó với mức độ này thêm phức tạp. Tôi bị cám dỗ chỉ phân tích là 2 trang web riêng biệt (thực tế là thời gian N_app không giống nhau ở mỗi trang web có thể khuyến khích phương pháp này). Tôi có thể sử dụng loại iii tổng của hình vuông ANOVA ở đây không?

Điều đó đã được đề xuất với tôi rằng phương pháp tạo mô hình hỗn hợp tuyến tính có thể là con đường phía trước nhưng tôi không quen sử dụng chúng.

Tôi sẽ hoan nghênh suy nghĩ của bạn về bất kỳ điều nào ở trên. Cảm ơn vì đã dành thời gian cho tôi.

Rory

Answer 1

Để trả lời câu hỏi đầu tiên của bạn về cách tốt nhất để thử nghiệm các giả định. Trong khi nỗ lực của bạn bằng cách sử dụng một thử nghiệm thống kê khác, được triển khai trong R, là hợp lý, tôi thực sự chỉ cần hình dung phân phối và xem liệu dữ liệu có đáp ứng các giả định ANOVA hay không. Cách tiếp cận này có vẻ hơi chủ quan, nhưng nó hoạt động trong hầu hết các trường hợp.

• độc lập, phân phối giống hệt (i.i.d) dữ liệu: đây là câu hỏi mà bạn có thể đã có câu trả lời dựa trên số lượng bạn biết về dữ liệu của mình. Có thể sử dụng phép thử chi-square để xác định tính độc lập (hay không).
• dữ liệu được phân phối bình thường: sử dụng biểu đồ lịch sử/QQ để kiểm tra. Dựa trên sự phân bố, tôi nghĩ rằng nó là hợp lý để sử dụng aov mặc dù phân phối hơi lưỡng cực.

(Dường như log-chuyển đổi giúp đỡ hơn nữa đáp ứng bình thường giả định. Đây là điều mà bạn có thể cân nhắc, nhất là đối với phân tích hạ lưu.)

par(mfrow=c(2,2)) 
plot(density(HF_NW_grass$Dry_tonnes_ha), col="red", main="Density") 
qqnorm(HF_NW_grass$Dry_tonnes_ha, col="red", main="qqplot") 
qqline(HF_NW_grass$Dry_tonnes_ha) 

DTH_trans <- log10(HF_NW_grass$Dry_tonnes_ha) 
plot(density(DTH_trans), col="blue", main="transformed density") 
qqnorm(DTH_trans, col="blue", main="transformed density") 
qqline(DTH_trans) 

Về câu hỏi thứ hai của bạn vào những gì là cách tốt nhất để xây dựng các biện pháp lặp lại vào mô hình là: Thật không may, rất khó để xác định mô hình "tốt nhất", nhưng dựa trên kiến ​​thức của tôi (chủ yếu thông qua dữ liệu lớn về gen), bạn có thể muốn sử dụng mô hình hiệu ứng hỗn hợp tuyến tính. Điều này có thể được thực hiện thông qua các gói lme4 R, ví dụ. Vì có vẻ như bạn đã biết cách xây dựng một mô hình tuyến tính trong R, bạn sẽ không gặp vấn đề gì với việc áp dụng các hàm lme4.

Câu hỏi thứ ba của bạn liên quan đến việc lồng hai biến là khó khăn. Nếu tôi là bạn, tôi sẽ bắt đầu với Site và Block như thể chúng là các yếu tố độc lập. Tuy nhiên, nếu bạn biết chúng không độc lập, có lẽ bạn nên lồng ghép chúng.

Tôi nghĩ câu hỏi và mối quan tâm của bạn khá cởi mở. Đề xuất của tôi là miễn là bạn có một lý giải chính đáng, hãy tiếp tục và tiếp tục.

Answer 2

Tôi đồng ý với @David C về việc sử dụng chẩn đoán bằng hình ảnh. lô QQ đơn giản nên làm việc

# dependent variable. 
par(mfrow=c(1,2)) 
qqnorm(dt[,dry_tonnes_ha]); qqline(dt[,dry_tonnes_ha], probs= c(0.15, 0.85)) 
qqnorm(log(dt[,dry_tonnes_ha])); qqline(log(dt[,dry_tonnes_ha]), probs= c(0.15, 0.85)) 

Việc chuyển đổi log trông hợp lý với tôi.Bạn cũng có thể thấy điều này từ ô mật độ, có chiều dài dài và có phần hai phương thức

par(mfrow=c(1,1)) 
plot(density(dt[,dry_tonnes_ha])) 

Bạn cũng có thể sử dụng ô xếp hàng (Buja et al, 2009) nếu muốn. Tôi không chắc họ cần thiết trong trường hợp này. Vignette provided

library(nullabor) 
# this may not be the best X variable. I'm not familiar with your data 
dt_l <- lineup(null_permute("dry_tonnes_ha"), dt) 
qplot(dry_tonnes_ha, treatment, data = dt_l) + facet_wrap(~ .sample) 

Đối với các giả định khác, bạn chỉ có thể sử dụng lô chẩn đoán tiêu chuẩn từ lm

lm2 <- lm(log(dry_tonnes_ha) ~ treatment * n_app + site/block, data = dt) 
plot(lm2) 

tôi không thấy bất cứ điều gì quá phiền hà trong các lô.