R中的2路随机块方差分析

Question

一直在尝试分析我的数据集，我相信我在正确的轨道上，但需要一些确认。我正在尝试分析同一条河流的几条河段的鱼类捕获率，并评估在 2 年的研究中使用的渔具类型的有效性。

我的数据包括：

1) 21 个不同的站点 2) 分类为“主动”或“被动”的采样技术 3) 2 年的数据收集，按月分开。

在研究过程中未对站点进行统一抽样。并非每个月都对它们进行采样，并非所有采样时间都相同，也不是采用相同的采样技术。我相信您可以将其归类为非重复测量，因为几乎没有两个采样周期是相同的。

我相信分析这些数据的正确方法是使用 2-way Randomized Block ANOVA。几个月将是分析中被阻止的内容。我得到了一些结果，但不确定使用的代码是否正确。

是否有人能够证明我使用的代码并确认/否认它确实是 R 中 2-Way Randomized Block 设计的正确代码？

Fish<-read.csv(file.choose(),header=TRUE)
Fish
FishLM<-lm(Caught.Hr ~ Site + Method + Site:Method,Fish)
anova(FishLM)

这里是一些示例数据：

Site     Month  Year    Device  Method  Hrs/Month   Caught  Caught/Hr  
Reach 01    5   2014    BS      Active  0.7            0    0  
Reach 01    6   2014    BS      Active  7.92           0    0  
Reach 01    7   2014    BS      Active  5.73           0    0  
Reach 01    8   2014    BS      Active  1.82           0    0  
Reach 01    9   2014    BS      Active  10.08          0    0  
Reach 01    10  2014    BS      Active  10.08          0    0  
Reach 01    11  2014    BS      Active  6.9            0    0  
Reach 02    3   2013    BS      Active  2.5            0    0  
Reach 02    4   2013    BS      Active  2.5            0    0  
Reach 02    5   2013    BS      Active  3.75           0    0  
Reach 02    6   2013    BS      Active  17.3           0    0  
Reach 02    7   2013    BS      Active  2.5            0    0  
Reach 02    8   2013    BS      Active  2.5            0    0  
Reach 02    9   2013    BS      Active  2.5            0    0  
Reach 02    10  2013    BS      Active  2.5            0    0  
Reach 02    11  2013    BS      Active  2.5            0    0  
Reach 03    3   2013    BS      Active  3              0    0  
Reach 03    4   2013    BS      Active  3              0    0  
Reach 03    5   2013    BS     Active   2.5            0    0  
Reach 03    6   2013    BS     Active   3.5            1    0.285714286  
Reach 03    7   2013    BS     Active   3              0    0  
Reach 03    8   2013    BS     Active   3              0    0  
Reach 03    9   2013    BS     Active   3              1    0.333333333  
Reach 03    10  2013    BS     Active   8.75           2    0.228571429  
Reach 03    11  2013    BS      Active  3              0    0  
Reach 04    3   2013    MT      Passive           
Reach 04    4   2013    MT      Passive           
Reach 04    5   2013    MT      Passive           
Reach 04    6   2013    MT      Passive 72             0    0  
Reach 04    7   2013    MT      Passive 120            2    0.016666667  
Reach 04    8   2013    MT      Passive 120            0    0  
Reach 04    9   2013    MT      Passive 72             0    0  
Reach 04    10  2013    MT      Passive           
Reach 04    11  2013    MT      Passive           
Reach 07    3   2014    MF      Passive           
Reach 07    4   2014    MF      Passive 96             7    0.072916667  
Reach 07    5   2014    MF      Passive 96             5    0.052083333  
Reach 07    6   2014    MF      Passive 96             8    0.083333333  
Reach 07    7   2014    MF      Passive 96             1    0.010416667  
Reach 07    8   2014    MF      Passive 96             1    0.010416667  
Reach 07    9   2014    MF      Passive 96             3    0.03125  
Reach 07    10  2014    MF      Passive 96            10    0.104166667  
Reach 07    11  2014    MF      Passive           

谢谢。

Answer 1

如果我正确阅读了您的问题（并且我不确定自己是否正确），您正在尝试检查您认为是正态分布的连续变量 Caught.Hr 的变化——因此是方差分析。另外，您有两个治疗效果：Site和Method，并且您每个月都有重复措施。

你的模型是这样的

$$Y_{ijk} = \mu + S_i + M_j + (SM){ij} + \epsilon{ijk}$$

在哪里
- Y_{ijk} 是在时间段 k 使用方法 j 在站点 i 的捕获率。
- mu 代表人口平均捕获率
- S_i表示每个站点的效果，
- M_j表示每种采样方式的效果，
- (SM)_{ij} 是交互效果，
- e_{ijk} 是随机变化

从您的描述中，我看不出您的阻碍因素是什么。听起来您的设计不平衡。根据您的描述，在我看来，您没有随机区组设计。您有一个包含两个因子的因子设计，这也是不平衡的。

但是是的，这会起作用：

FishLM<-lm(Caught.Hr ~ Site + Method + Site:Method,Fish)
anova(FishLM)

编辑：

根据您的数据，我认为我上面所说的内容是有效的。尽管我确实担心您使用的是方差分析。这似乎是计数数据，即泊松不是正态分布。例如：

# This has problems based on hours / obs
Fishglm <- glm(Caught ~ Site + Method + Site:Method, data=Fish, 
  family= poisson(link = "log"))
# could use neg-binomial on the rate instead.
library(MASS)
Fishnb <- glm.nb(Caught.Hr ~ Site + Method + Site:Method, data=Fish)