Cumulative Precipitation
Both sites received less precipitation in 2017 than 2016.
Animal Gains
HREC Animal Gains
HrecCalves <- HrecCalves %>% mutate(Pasture=gsub(" ","",Pasture))
## Warning: package 'bindrcpp' was built under R version 3.4.3
HrecCows <- HrecCows %>% mutate(Pasture=gsub(" ","",Pasture))
HrecSheep <- HrecSheep %>% mutate(Pasture=gsub(" ","",Pasture))
ADG_HREC_Calves <- aov(avgDday~Year + Error(Pasture), data=HrecCalves)
summary(ADG_HREC_Calves) # f value 4.466 Pr(>F) 0.0363
##
## Error: Pasture
## Df Sum Sq Mean Sq F value Pr(>F)
## Year 1 0.6997 0.6997 1.063 0.49
## Residuals 1 0.6580 0.6580
##
## Error: Within
## Df Sum Sq Mean Sq F value Pr(>F)
## Year 1 0.38 0.3804 4.466 0.0363 *
## Residuals 143 12.18 0.0852
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
ADG_HREC_Cows <- aov(avgDday~Year + Error(Pasture), data =HrecCows)
summary(ADG_HREC_Cows) # f value 3.827 Pr(>F) 0.0523
##
## Error: Pasture
## Df Sum Sq Mean Sq F value Pr(>F)
## Year 1 1.0685 1.0685 3.615 0.308
## Residuals 1 0.2956 0.2956
##
## Error: Within
## Df Sum Sq Mean Sq F value Pr(>F)
## Year 1 1.43 1.4334 3.827 0.0523 .
## Residuals 150 56.19 0.3746
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
ADG_Hrec_Sheep <-aov(avgDday~Year + Error(Pasture), data=HrecSheep)
summary(ADG_Hrec_Sheep) # f value 29.95 Pr(>F) 5.51e-08
##
## Error: Pasture
## Df Sum Sq Mean Sq F value Pr(>F)
## Year 1 0.05488 0.05488 0.199 0.733
## Residuals 1 0.27568 0.27568
##
## Error: Within
## Df Sum Sq Mean Sq F value Pr(>F)
## Year 1 0.274 0.27376 29.95 5.51e-08 ***
## Residuals 1090 9.965 0.00914
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
HrecCalves <- HrecCalves %>% mutate(Pasture=gsub(" ","",Pasture))
HrecCalves$Year <- factor(HrecCalves$Year)
HrecCows <- HrecCows %>% mutate(Pasture=gsub(" ","",Pasture))
HrecCows$Year <- factor(HrecCows$Year)
HrecSheep <- HrecSheep %>% mutate(Pasture=gsub(" ","",Pasture))
HrecSheep$Year <- factor(HrecSheep$Year)
ADG_HREC_Calves2 <- lmer(avgDday~Year + 0 + (1|Pasture), data=HrecCalves)
Mult_HREC_Calves2 <- glht(ADG_HREC_Calves2, linfct=c("Year2016 = 0",
"Year2017 = 0"))
summary(Mult_HREC_Calves2)
##
## Simultaneous Tests for General Linear Hypotheses
##
## Fit: lmer(formula = avgDday ~ Year + 0 + (1 | Pasture), data = HrecCalves)
##
## Linear Hypotheses:
## Estimate Std. Error z value Pr(>|z|)
## Year2016 == 0 2.99256 0.07031 42.56 <1e-10 ***
## Year2017 == 0 3.09562 0.07014 44.13 <1e-10 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## (Adjusted p values reported -- single-step method)
confint(ADG_HREC_Calves2)
## Computing profile confidence intervals ...
## 2.5 % 97.5 %
## .sig01 0.02663542 0.2721540
## .sigma 0.26022606 0.3279507
## Year2016 2.83880840 3.1479997
## Year2017 2.94142296 3.2502626
ADG_HREC_Cows2 <- lmer(avgDday~Year + 0 + (1|Pasture), data=HrecCows)
Mult_HREC_Cows2 <- glht(ADG_HREC_Cows2, linfct=c("Year2016 = 0",
"Year2017 = 0"))
summary(Mult_HREC_Cows2)
##
## Simultaneous Tests for General Linear Hypotheses
##
## Fit: lmer(formula = avgDday ~ Year + 0 + (1 | Pasture), data = HrecCows)
##
## Linear Hypotheses:
## Estimate Std. Error z value Pr(>|z|)
## Year2016 == 0 0.64285 0.08488 7.574 < 1e-10 ***
## Year2017 == 0 0.45370 0.08556 5.303 2.28e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## (Adjusted p values reported -- single-step method)
confint(ADG_HREC_Cows2)
## Computing profile confidence intervals ...
## Warning in optwrap(optimizer, par = thopt, fn = mkdevfun(rho, 0L), lower
## = fitted@lower): convergence code 3 from bobyqa: bobyqa -- a trust region
## step failed to reduce q
## 2.5 % 97.5 %
## .sig01 0.0000000 0.2764354
## .sigma 0.5470750 0.6856556
## Year2016 0.4706808 0.8205275
## Year2017 0.2786964 0.6300638
ADG_HREC_Sheep2 <- lmer(avgDday~Year + 0 + (1|Pasture), data=HrecSheep)
Mult_HREC_Sheep2 <- glht(ADG_HREC_Sheep2, linfct=c("Year2016 = 0",
"Year2017 = 0"))
summary(Mult_HREC_Sheep2)
##
## Simultaneous Tests for General Linear Hypotheses
##
## Fit: lmer(formula = avgDday ~ Year + 0 + (1 | Pasture), data = HrecSheep)
##
## Linear Hypotheses:
## Estimate Std. Error z value Pr(>|z|)
## Year2016 == 0 0.25761 0.01255 20.52 <1e-10 ***
## Year2017 == 0 0.28930 0.01263 22.90 <1e-10 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## (Adjusted p values reported -- single-step method)
confint(ADG_HREC_Sheep2)
## Computing profile confidence intervals ...
## 2.5 % 97.5 %
## .sig01 0.007734052 0.05092279
## .sigma 0.091696304 0.09972491
## Year2016 0.229388738 0.28582508
## Year2017 0.261004367 0.31759827
Calves and sheep are statistically different. Calves performed better in 2017. Sheep are slightly better in 2017? There is very little variation combined with more sheep which seems to account for the significant difference.
Still need to update graphs!!!
HREC ADG Violins
HREC ADG +/- sd
ggplot(data=subset(REC_ADG, Location=="HREC"), aes(x=Type, y=mean, color=Type)) +
geom_errorbar(aes(ymin=mean-sd, ymax=mean+sd), width =0.25, lwd=1.5)+
theme_bw() +
geom_hline(yintercept = 0, linetype="dashed") +
geom_point(size=4) +
labs(y="Average Daily Gain (lbs/day)") +
facet_wrap(~Year, ncol = 2) +
theme(axis.text.x = element_text(angle=33, hjust=1))
HREC ADG +/- se
ggplot(data=subset(REC_ADG, Location=="HREC"), aes(x=Type, y=mean, color=Type)) +
geom_errorbar(aes(ymin=mean-se, ymax=mean+se), width =0.25, lwd=1.5)+
theme_bw() +
geom_hline(yintercept = 0, linetype="dashed") +
geom_point(size=4) +
labs(y="Average Daily Gain (lbs/day)") +
facet_wrap(~Year, ncol = 2) +
theme(axis.text.x = element_text(angle=33, hjust=1))
CGREC Animal Gains
#comparing treatments against each other
ADG_CGREC_Calves <- lmer(avgDday~Treatment + (1|Pasture), data=CgrecCalves)
Mult_CREC_Calves <- glht(ADG_CGREC_Calves, linfct=mcp(Treatment = "Tukey"))
summary(Mult_CREC_Calves)
##
## Simultaneous Tests for General Linear Hypotheses
##
## Multiple Comparisons of Means: Tukey Contrasts
##
##
## Fit: lmer(formula = avgDday ~ Treatment + (1 | Pasture), data = CgrecCalves)
##
## Linear Hypotheses:
## Estimate Std. Error z value Pr(>|z|)
## SpringOnly - NoFire == 0 0.09410 0.05204 1.808 0.167
## SpringSummer - NoFire == 0 0.06183 0.05074 1.218 0.442
## SpringSummer - SpringOnly == 0 -0.03227 0.05109 -0.632 0.803
## (Adjusted p values reported -- single-step method)
ADG_CGREC_Cows <- lmer(avgDday~Treatment + (1|Pasture), data=CgrecCows)
Mult_CGREC_Cows <- glht(ADG_CGREC_Cows, linfct=mcp(Treatment = "Tukey"))
summary(Mult_CGREC_Cows) #both fire treatments are diffent than the no fire
##
## Simultaneous Tests for General Linear Hypotheses
##
## Multiple Comparisons of Means: Tukey Contrasts
##
##
## Fit: lmer(formula = avgDday ~ Treatment + (1 | Pasture), data = CgrecCows)
##
## Linear Hypotheses:
## Estimate Std. Error z value Pr(>|z|)
## SpringOnly - NoFire == 0 1.18656 0.46610 2.546 0.0293 *
## SpringSummer - NoFire == 0 1.23397 0.46566 2.650 0.0220 *
## SpringSummer - SpringOnly == 0 0.04741 0.46587 0.102 0.9943
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## (Adjusted p values reported -- single-step method)
#get the adg and confidence interval for each treatment
ADG_CGREC_Calves2 <- lmer(avgDday~Treatment + 0 + (1|Pasture), data=CgrecCalves)
Mult_CREC_Calves2 <- glht(ADG_CGREC_Calves2, linfct=c("TreatmentNoFire = 0",
"TreatmentSpringSummer = 0",
"TreatmentSpringOnly = 0"))
summary(Mult_CREC_Calves2) #different than zero
##
## Simultaneous Tests for General Linear Hypotheses
##
## Fit: lmer(formula = avgDday ~ Treatment + 0 + (1 | Pasture), data = CgrecCalves)
##
## Linear Hypotheses:
## Estimate Std. Error z value Pr(>|z|)
## TreatmentNoFire == 0 2.58929 0.03655 70.83 <2e-16 ***
## TreatmentSpringSummer == 0 2.65112 0.03520 75.32 <2e-16 ***
## TreatmentSpringOnly == 0 2.68339 0.03703 72.46 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## (Adjusted p values reported -- single-step method)
confint(ADG_CGREC_Calves2)
## Computing profile confidence intervals ...
## Warning in nextpar(mat, cc, i, delta, lowcut, upcut): unexpected decrease
## in profile: using minstep
## Warning in FUN(X[[i]], ...): non-monotonic profile for .sig01
## Warning in confint.thpr(pp, level = level, zeta = zeta): bad spline fit
## for .sig01: falling back to linear interpolation
## 2.5 % 97.5 %
## .sig01 0.0000000 0.0781369
## .sigma 0.2790288 0.3285863
## TreatmentNoFire 2.5240188 2.6566023
## TreatmentSpringOnly 2.6190493 2.7443049
## TreatmentSpringSummer 2.5941232 2.7100873
fixef(ADG_CGREC_Calves2)
## TreatmentNoFire TreatmentSpringOnly TreatmentSpringSummer
## 2.589287 2.683385 2.651116
ADG_CGREC_Cows2 <- lmer(avgDday~Treatment + 0 + (1|Pasture), data=CgrecCows)
Mult_CGREC_Cows2 <- glht(ADG_CGREC_Cows2, linfct=c("TreatmentNoFire = 0",
"TreatmentSpringSummer = 0",
"TreatmentSpringOnly = 0"))
summary(Mult_CGREC_Cows2)
##
## Simultaneous Tests for General Linear Hypotheses
##
## Fit: lmer(formula = avgDday ~ Treatment + 0 + (1 | Pasture), data = CgrecCows)
##
## Linear Hypotheses:
## Estimate Std. Error z value Pr(>|z|)
## TreatmentNoFire == 0 -0.5045 0.3294 -1.531 0.3316
## TreatmentSpringSummer == 0 0.7294 0.3291 2.216 0.0779 .
## TreatmentSpringOnly == 0 0.6820 0.3297 2.068 0.1114
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## (Adjusted p values reported -- single-step method)
confint(ADG_CGREC_Cows2)
## Computing profile confidence intervals ...
## 2.5 % 97.5 %
## .sig01 0.38818804 0.8997925
## .sigma 0.42223219 0.5003632
## TreatmentNoFire -1.11142191 0.1024986
## TreatmentSpringOnly 0.07450877 1.2894870
## TreatmentSpringSummer 0.12298592 1.3357239
fixef(ADG_CGREC_Cows2)
## TreatmentNoFire TreatmentSpringOnly TreatmentSpringSummer
## -0.5045274 0.6820291 0.7294417
No difference in the calves, but both burn treatments are different than the no burn treatment.
CGREC ADG Violins
ggplot(data=subset(Livestock, Location=="CGREC"), aes(x=Type, y=avgDday)) + theme_bw(16) +
geom_hline(yintercept = 0) +
geom_violin() +
geom_jitter(aes(color=Pasture), width=0.3,
alpha=0.75, show.legend=FALSE) +
labs(y="Average Daily Gain (lbs/day)") +
facet_wrap(~Treatment, ncol = 3) +
theme(axis.text.x = element_text(angle=33, hjust=1))
CGREC ADG +/- SD
ggplot(data=subset(REC_ADG, Location=="CGREC"), aes(x=Type, y=mean, color=Type)) +
geom_errorbar(aes(ymin=mean-sd, ymax=mean+sd), width =0.25, lwd=1.5)+
theme_bw() +
geom_hline(yintercept = 0, linetype="dashed") +
geom_point(size=4) +
labs(y="Average Daily Gain (lbs/day)") +
facet_wrap(~Treatment, ncol = 3) +
theme(axis.text.x = element_text(angle=33, hjust=1))
CGREC ADG +/- SE
ggplot(data=subset(REC_ADG, Location=="CGREC"), aes(x=Type, y=mean, color=Type)) +
geom_errorbar(aes(ymin=mean-se, ymax=mean+se), width =0.25, lwd=1.5)+
theme_bw() +
geom_hline(yintercept = 0, linetype="dashed") +
geom_point(size=4) +
labs(y="Average Daily Gain (lbs/day)") +
facet_wrap(~Treatment, ncol = 3) +
theme(axis.text.x = element_text(angle=33, hjust=1))
Confidence intervals for Fecal and Protein
I think the bootstrap confidence interval strategy is not working for the fecal counts and protein due to a lack of data points that go into them. I tried using each month as an individual data point to increase the sample size, but that did not do much. The monthly mean difference is how I first tried to test this and incorporate the repeated monthly sampling.
Fecal.HREC.null <- lmer(Fecal ~ 1 + (1|Pasture/Month),
data=subset(ProGO, REC=="HREC"), REML=FALSE)
Fecal.HREC.p <- lmer(Fecal ~ Patch + (1|Pasture/Month),
data=subset(ProGO, REC=="HREC"), REML=FALSE)
Fecal.HREC.pt <- lmer(Fecal ~ Patch + Treatment + (1|Pasture/Month),
data=subset(ProGO, REC=="HREC"), REML=FALSE)
anova(Fecal.HREC.null, Fecal.HREC.p, Fecal.HREC.pt) %>%
pander(caption="comparison of the three models")
comparison of the three models (continued below)
| Fecal.HREC.null |
4 |
1975 |
1990 |
-983.3 |
1967 |
NA |
NA |
| Fecal.HREC.p |
5 |
1934 |
1954 |
-962.1 |
1924 |
42.44 |
1 |
| Fecal.HREC.pt |
6 |
1927 |
1950 |
-957.3 |
1915 |
9.628 |
1 |
| Fecal.HREC.null |
NA |
| Fecal.HREC.p |
7.292e-11 |
| Fecal.HREC.pt |
0.001916 |
anova(Fecal.HREC.null, Fecal.HREC.pt) %>%
pander(caption="null vs model with patch and treatment")
null vs model with patch and treatment (continued below)
| Fecal.HREC.null |
4 |
1975 |
1990 |
-983.3 |
1967 |
NA |
NA |
| Fecal.HREC.pt |
6 |
1927 |
1950 |
-957.3 |
1915 |
52.07 |
2 |
| Fecal.HREC.null |
NA |
| Fecal.HREC.pt |
4.941e-12 |
anova(Fecal.HREC.null, Fecal.HREC.p) %>%
pander(caption="null vs model with patch")
null vs model with patch (continued below)
| Fecal.HREC.null |
4 |
1975 |
1990 |
-983.3 |
1967 |
NA |
NA |
| Fecal.HREC.p |
5 |
1934 |
1954 |
-962.1 |
1924 |
42.44 |
1 |
| Fecal.HREC.null |
NA |
| Fecal.HREC.p |
7.292e-11 |
Fecal.HREC.INT <- lmer(Fecal ~ Patch * Treatment + (1|Pasture/Month),
data=subset(ProGO, REC=="HREC"), REML=FALSE)
anova(Fecal.HREC.pt, Fecal.HREC.INT)
## Data: subset(ProGO, REC == "HREC")
## Models:
## Fecal.HREC.pt: Fecal ~ Patch + Treatment + (1 | Pasture/Month)
## Fecal.HREC.INT: Fecal ~ Patch * Treatment + (1 | Pasture/Month)
## Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq)
## Fecal.HREC.pt 6 1926.6 1950.3 -957.29 1914.6
## Fecal.HREC.INT 7 1925.8 1953.6 -955.92 1911.8 2.7382 1 0.09798
##
## Fecal.HREC.pt
## Fecal.HREC.INT .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary(Fecal.HREC.pt)
## Linear mixed model fit by maximum likelihood ['lmerMod']
## Formula: Fecal ~ Patch + Treatment + (1 | Pasture/Month)
## Data: subset(ProGO, REC == "HREC")
##
## AIC BIC logLik deviance df.resid
## 1926.6 1950.3 -957.3 1914.6 381
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.8643 -0.6093 -0.1039 0.4130 5.4006
##
## Random effects:
## Groups Name Variance Std.Dev.
## Month:Pasture (Intercept) 9.807e-01 9.903e-01
## Pasture (Intercept) 1.280e-13 3.578e-07
## Residual 7.755e+00 2.785e+00
## Number of obs: 387, groups: Month:Pasture, 18; Pasture, 6
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 4.0490 0.4590 8.821
## PatchU -2.2314 0.3307 -6.748
## TreatmentSheep 2.2983 0.5463 4.207
##
## Correlation of Fixed Effects:
## (Intr) PatchU
## PatchU -0.544
## TreatmntShp -0.591 0.000
confint(Fecal.HREC.pt)
## Computing profile confidence intervals ...
## Warning in zeta(shiftpar, start = opt[seqpar1][-w]): slightly lower
## deviances (diff=-2.27374e-13) detected
## Warning in nextpar(mat, cc, i, delta, lowcut, upcut): Last two rows have
## identical or NA .zeta values: using minstep
## Warning in nextpar(mat, cc, i, delta, lowcut, upcut): Last two rows have
## identical or NA .zeta values: using minstep
## Warning in zetafun(np, ns): slightly lower deviances (diff=-2.27374e-13)
## detected
## Warning in FUN(X[[i]], ...): non-monotonic profile for .sig02
## Warning in confint.thpr(pp, level = level, zeta = zeta): bad spline fit
## for .sig02: falling back to linear interpolation
## 2.5 % 97.5 %
## .sig01 0.6117623 1.5615556
## .sig02 0.0000000 0.9250144
## .sigma 2.5953345 2.9985164
## (Intercept) 3.1236080 4.9747728
## PatchU -2.8814333 -1.5816383
## TreatmentSheep 1.1584913 3.4335414
Fecal.CGREC.null <- lmer(Fecal ~ 1 + (1|Pasture/Month),
data=subset(ProGO, REC=="CGREC"), REML=FALSE)
Fecal.CGREC.p <- lmer(Fecal ~ Patch + (1|Pasture/Month),
data=subset(ProGO, REC=="CGREC"), REML=FALSE)
Fecal.CGREC.pt <- lmer(Fecal ~ Patch + Treatment + (1|Pasture/Month),
data=subset(ProGO, REC=="CGREC"), REML=FALSE)
anova(Fecal.CGREC.null, Fecal.CGREC.p, Fecal.CGREC.pt) %>%
pander(caption="comparison of the three models")
comparison of the three models (continued below)
| Fecal.CGREC.null |
4 |
2765 |
2784 |
-1379 |
2757 |
NA |
NA |
| Fecal.CGREC.p |
5 |
2589 |
2613 |
-1290 |
2579 |
177.7 |
1 |
| Fecal.CGREC.pt |
6 |
2587 |
2614 |
-1287 |
2575 |
4.955 |
1 |
| Fecal.CGREC.null |
NA |
| Fecal.CGREC.p |
1.539e-40 |
| Fecal.CGREC.pt |
0.02602 |
anova(Fecal.CGREC.null,Fecal.CGREC.pt) %>%
pander(caption="null vs patch and treatment")
null vs patch and treatment (continued below)
| Fecal.CGREC.null |
4 |
2765 |
2784 |
-1379 |
2757 |
NA |
NA |
| Fecal.CGREC.pt |
6 |
2587 |
2614 |
-1287 |
2575 |
182.7 |
2 |
| Fecal.CGREC.null |
NA |
| Fecal.CGREC.pt |
2.171e-40 |
anova(Fecal.CGREC.null, Fecal.CGREC.p) %>%
pander(caption="null vs patch only")
null vs patch only (continued below)
| Fecal.CGREC.null |
4 |
2765 |
2784 |
-1379 |
2757 |
NA |
NA |
| Fecal.CGREC.p |
5 |
2589 |
2613 |
-1290 |
2579 |
177.7 |
1 |
| Fecal.CGREC.null |
NA |
| Fecal.CGREC.p |
1.539e-40 |
anova(Fecal.CGREC.pt, Fecal.CGREC.p) %>%
pander(caption="patch and treatment vs patch only")
patch and treatment vs patch only (continued below)
| Fecal.CGREC.p |
5 |
2589 |
2613 |
-1290 |
2579 |
NA |
NA |
| Fecal.CGREC.pt |
6 |
2587 |
2614 |
-1287 |
2575 |
4.955 |
1 |
| Fecal.CGREC.p |
NA |
| Fecal.CGREC.pt |
0.02602 |
Fecal.CGREC.INT <- lmer(Fecal ~ Patch * Treatment + (1|Pasture/Month),
data=subset(ProGO, REC=="CGREC"), REML=FALSE)
anova(Fecal.CGREC.pt, Fecal.CGREC.INT)
## Data: subset(ProGO, REC == "CGREC")
## Models:
## Fecal.CGREC.pt: Fecal ~ Patch + Treatment + (1 | Pasture/Month)
## Fecal.CGREC.INT: Fecal ~ Patch * Treatment + (1 | Pasture/Month)
## Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq)
## Fecal.CGREC.pt 6 2586.5 2614.3 -1287.3 2574.5
## Fecal.CGREC.INT 7 2588.5 2620.9 -1287.3 2574.5 0.0021 1 0.9638
summary(Fecal.CGREC.pt)
## Linear mixed model fit by maximum likelihood ['lmerMod']
## Formula: Fecal ~ Patch + Treatment + (1 | Pasture/Month)
## Data: subset(ProGO, REC == "CGREC")
##
## AIC BIC logLik deviance df.resid
## 2586.5 2614.3 -1287.3 2574.5 751
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.2632 -0.5794 -0.2970 0.3760 7.7836
##
## Random effects:
## Groups Name Variance Std.Dev.
## Month:Pasture (Intercept) 0.07048 0.2655
## Pasture (Intercept) 0.00000 0.0000
## Residual 1.71050 1.3079
## Number of obs: 757, groups: Month:Pasture, 24; Pasture, 8
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 2.6878 0.1515 17.737
## PatchU -1.8529 0.1295 -14.311
## TreatmentSpring Only -0.3385 0.1448 -2.337
##
## Correlation of Fixed Effects:
## (Intr) PatchU
## PatchU -0.743
## TrtmntSprnO -0.515 0.062
confint(Fecal.CGREC.pt)
## Computing profile confidence intervals ...
## 2.5 % 97.5 %
## .sig01 0.1380031 0.42457752
## .sig02 0.0000000 0.19460993
## .sigma 1.2436480 1.37779835
## (Intercept) 2.3886292 2.98842573
## PatchU -2.1069856 -1.59878928
## TreatmentSpring Only -0.6353726 -0.04389973
summary(Fecal.CGREC.p)
## Linear mixed model fit by maximum likelihood ['lmerMod']
## Formula: Fecal ~ Patch + (1 | Pasture/Month)
## Data: subset(ProGO, REC == "CGREC")
##
## AIC BIC logLik deviance df.resid
## 2589.5 2612.6 -1289.7 2579.5 752
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.3702 -0.5547 -0.3021 0.3328 7.8266
##
## Random effects:
## Groups Name Variance Std.Dev.
## Month:Pasture (Intercept) 9.851e-02 3.139e-01
## Pasture (Intercept) 1.996e-47 4.468e-24
## Residual 1.711e+00 1.308e+00
## Number of obs: 757, groups: Month:Pasture, 24; Pasture, 8
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 2.5087 0.1344 18.67
## PatchU -1.8380 0.1293 -14.21
##
## Correlation of Fixed Effects:
## (Intr)
## PatchU -0.804
confint(Fecal.CGREC.p)
## Computing profile confidence intervals ...
## Warning in zeta(shiftpar, start = opt[seqpar1][-w]): slightly lower
## deviances (diff=-4.54747e-13) detected
## Warning in nextpar(mat, cc, i, delta, lowcut, upcut): Last two rows have
## identical or NA .zeta values: using minstep
## Warning in nextpar(mat, cc, i, delta, lowcut, upcut): Last two rows have
## identical or NA .zeta values: using minstep
## Warning in zetafun(np, ns): slightly lower deviances (diff=-4.54747e-13)
## detected
## Warning in FUN(X[[i]], ...): non-monotonic profile for .sig02
## Warning in confint.thpr(pp, level = level, zeta = zeta): bad spline fit
## for .sig02: falling back to linear interpolation
## 2.5 % 97.5 %
## .sig01 0.1818811 0.4831699
## .sig02 0.0000000 0.2972714
## .sigma 1.2437317 1.3779110
## (Intercept) 2.2444090 2.7746637
## PatchU -2.0919826 -1.5838729
Protein
Pro.HREC.null <- lmer(PROTEIN ~ 1 + (1|Pasture/Month),
data=subset(ProGO, REC=="HREC"), REML=FALSE)
Pro.HREC.p <- lmer(PROTEIN ~ Patch + (1|Pasture/Month),
data=subset(ProGO, REC=="HREC"), REML=FALSE)
Pro.HREC.pt <- lmer(PROTEIN ~ Patch + Treatment + (1|Pasture/Month),
data=subset(ProGO, REC=="HREC"), REML=FALSE)
anova(Pro.HREC.null, Pro.HREC.p, Pro.HREC.pt) %>%
pander(caption="comparison of the three models")
comparison of the three models (continued below)
| Pro.HREC.null |
4 |
1708 |
1724 |
-850.1 |
1700 |
NA |
NA |
| Pro.HREC.p |
5 |
1707 |
1727 |
-848.5 |
1697 |
3.169 |
1 |
| Pro.HREC.pt |
6 |
1703 |
1726 |
-845.3 |
1691 |
6.373 |
1 |
| Pro.HREC.null |
NA |
| Pro.HREC.p |
0.07505 |
| Pro.HREC.pt |
0.01158 |
#Not different than the null model
anova(Pro.HREC.null, Pro.HREC.p) %>%
pander(caption="null vs patch only")
null vs patch only (continued below)
| Pro.HREC.null |
4 |
1708 |
1724 |
-850.1 |
1700 |
NA |
NA |
| Pro.HREC.p |
5 |
1707 |
1727 |
-848.5 |
1697 |
3.169 |
1 |
| Pro.HREC.null |
NA |
| Pro.HREC.p |
0.07505 |
anova(Pro.HREC.null, Pro.HREC.pt) %>%
pander(caption="null vs patch and treatment")
null vs patch and treatment (continued below)
| Pro.HREC.null |
4 |
1708 |
1724 |
-850.1 |
1700 |
NA |
NA |
| Pro.HREC.pt |
6 |
1703 |
1726 |
-845.3 |
1691 |
9.542 |
2 |
| Pro.HREC.null |
NA |
| Pro.HREC.pt |
0.00847 |
anova(Pro.HREC.p, Pro.HREC.pt) %>%
pander(caption="patch vs patch and treatment")
patch vs patch and treatment
| Pro.HREC.p |
5 |
1707 |
1727 |
-848.5 |
1697 |
NA |
NA |
NA |
| Pro.HREC.pt |
6 |
1703 |
1726 |
-845.3 |
1691 |
6.373 |
1 |
0.01158 |
Pro.HREC.INT <- lmer(PROTEIN ~ Patch * Treatment + (1|Pasture/Month),
data=subset(ProGO, REC=="HREC"), REML=FALSE)
anova(Pro.HREC.pt, Pro.HREC.INT)
## Data: subset(ProGO, REC == "HREC")
## Models:
## Pro.HREC.pt: PROTEIN ~ Patch + Treatment + (1 | Pasture/Month)
## Pro.HREC.INT: PROTEIN ~ Patch * Treatment + (1 | Pasture/Month)
## Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq)
## Pro.HREC.pt 6 1702.6 1726.3 -845.28 1690.6
## Pro.HREC.INT 7 1704.4 1732.1 -845.19 1690.4 0.1761 1 0.6748
summary(Pro.HREC.pt)
## Linear mixed model fit by maximum likelihood ['lmerMod']
## Formula: PROTEIN ~ Patch + Treatment + (1 | Pasture/Month)
## Data: subset(ProGO, REC == "HREC")
##
## AIC BIC logLik deviance df.resid
## 1702.6 1726.3 -845.3 1690.6 381
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.4975 -0.6211 -0.1070 0.4765 5.3782
##
## Random effects:
## Groups Name Variance Std.Dev.
## Month:Pasture (Intercept) 0.26903 0.5187
## Pasture (Intercept) 0.07008 0.2647
## Residual 4.41918 2.1022
## Number of obs: 387, groups: Month:Pasture, 18; Pasture, 6
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 10.0808 0.3335 30.226
## PatchU -0.4386 0.2496 -1.757
## TreatmentSheep -1.3074 0.3904 -3.349
##
## Correlation of Fixed Effects:
## (Intr) PatchU
## PatchU -0.566
## TreatmntShp -0.581 0.000
confint(Pro.HREC.pt) #no protein difference at patch, sheep lower than cattle
## Computing profile confidence intervals ...
## 2.5 % 97.5 %
## .sig01 0.1687074 0.94355440
## .sig02 0.0000000 0.87170367
## .sigma 1.9592320 2.26340034
## (Intercept) 9.3657468 10.79909521
## PatchU -0.9291927 0.05209052
## TreatmentSheep -2.2147505 -0.40692034
fixef(Pro.HREC.pt)
## (Intercept) PatchU TreatmentSheep
## 10.080797 -0.438619 -1.307373
summary(Pro.HREC.p)
## Linear mixed model fit by maximum likelihood ['lmerMod']
## Formula: PROTEIN ~ Patch + (1 | Pasture/Month)
## Data: subset(ProGO, REC == "HREC")
##
## AIC BIC logLik deviance df.resid
## 1706.9 1726.7 -848.5 1696.9 382
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.4440 -0.6055 -0.1066 0.4862 5.4384
##
## Random effects:
## Groups Name Variance Std.Dev.
## Month:Pasture (Intercept) 0.2690 0.5186
## Pasture (Intercept) 0.4964 0.7046
## Residual 4.4199 2.1024
## Number of obs: 387, groups: Month:Pasture, 18; Pasture, 6
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 9.4290 0.3805 24.782
## PatchU -0.4456 0.2497 -1.784
##
## Correlation of Fixed Effects:
## (Intr)
## PatchU -0.496
confint(Pro.HREC.p)
## Computing profile confidence intervals ...
## 2.5 % 97.5 %
## .sig01 0.1671093 0.99173263
## .sig02 0.1691537 1.57480330
## .sigma 1.9593667 2.26360723
## (Intercept) 8.5959614 10.26074984
## PatchU -0.9363244 0.04522755
Pro.CGREC.null <- lmer(PROTEIN ~ 1 + (1|Pasture/Month),
data=subset(ProGO, REC=="CGREC"), REML=FALSE)
Pro.CGREC.p <- lmer(PROTEIN ~ Patch + (1|Pasture/Month),
data=subset(ProGO, REC=="CGREC"), REML=FALSE)
Pro.CGREC.pt <- lmer(PROTEIN ~ Patch + Treatment + (1|Pasture/Month),
data=subset(ProGO, REC=="CGREC"), REML=FALSE)
anova(Pro.CGREC.null, Pro.CGREC.p, Pro.CGREC.pt) %>%
pander(caption="comparison of the three models")
comparison of the three models (continued below)
| Pro.CGREC.null |
4 |
3383 |
3402 |
-1688 |
3375 |
NA |
NA |
| Pro.CGREC.p |
5 |
3160 |
3183 |
-1575 |
3150 |
225.1 |
1 |
| Pro.CGREC.pt |
6 |
3162 |
3190 |
-1575 |
3150 |
0.0004245 |
1 |
| Pro.CGREC.null |
NA |
| Pro.CGREC.p |
6.906e-51 |
| Pro.CGREC.pt |
0.9836 |
#.P IS THE BETTER MODEL
anova(Pro.CGREC.null, Pro.CGREC.p) %>%
pander(caption=" null vs patch")
null vs patch (continued below)
| Pro.CGREC.null |
4 |
3383 |
3402 |
-1688 |
3375 |
NA |
NA |
| Pro.CGREC.p |
5 |
3160 |
3183 |
-1575 |
3150 |
225.1 |
1 |
| Pro.CGREC.null |
NA |
| Pro.CGREC.p |
6.906e-51 |
anova(Pro.CGREC.null, Pro.CGREC.pt) %>%
pander(caption="null vs patch and treatment")
null vs patch and treatment (continued below)
| Pro.CGREC.null |
4 |
3383 |
3402 |
-1688 |
3375 |
NA |
NA |
| Pro.CGREC.pt |
6 |
3162 |
3190 |
-1575 |
3150 |
225.1 |
2 |
| Pro.CGREC.null |
NA |
| Pro.CGREC.pt |
1.304e-49 |
anova(Pro.CGREC.p, Pro.CGREC.pt) %>%
pander(caption="patch vs patch and treatment")
patch vs patch and treatment (continued below)
| Pro.CGREC.p |
5 |
3160 |
3183 |
-1575 |
3150 |
NA |
NA |
| Pro.CGREC.pt |
6 |
3162 |
3190 |
-1575 |
3150 |
0.0004245 |
1 |
| Pro.CGREC.p |
NA |
| Pro.CGREC.pt |
0.9836 |
Pro.CGREC.INT <- lmer(PROTEIN ~ Patch * Treatment + (1|Pasture/Month),
data=subset(ProGO, REC=="CGREC"), REML=FALSE)
anova(Pro.CGREC.pt, Pro.CGREC.INT)
## Data: subset(ProGO, REC == "CGREC")
## Models:
## Pro.CGREC.pt: PROTEIN ~ Patch + Treatment + (1 | Pasture/Month)
## Pro.CGREC.INT: PROTEIN ~ Patch * Treatment + (1 | Pasture/Month)
## Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq)
## Pro.CGREC.pt 6 3162.0 3189.7 -1575.0 3150.0
## Pro.CGREC.INT 7 3137.5 3169.9 -1561.8 3123.5 26.428 1 2.735e-07
##
## Pro.CGREC.pt
## Pro.CGREC.INT ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary(Pro.CGREC.pt)
## Linear mixed model fit by maximum likelihood ['lmerMod']
## Formula: PROTEIN ~ Patch + Treatment + (1 | Pasture/Month)
## Data: subset(ProGO, REC == "CGREC")
##
## AIC BIC logLik deviance df.resid
## 3162.0 3189.7 -1575.0 3150.0 751
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.1770 -0.6828 -0.1490 0.5190 4.9338
##
## Random effects:
## Groups Name Variance Std.Dev.
## Month:Pasture (Intercept) 0.377739 0.6146
## Pasture (Intercept) 0.008836 0.0940
## Residual 3.583823 1.8931
## Number of obs: 757, groups: Month:Pasture, 24; Pasture, 8
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 12.188125 0.263938 46.178
## PatchU -3.040369 0.187688 -16.199
## TreatmentSpring Only 0.006072 0.294560 0.021
##
## Correlation of Fixed Effects:
## (Intr) PatchU
## PatchU -0.619
## TrtmntSprnO -0.581 0.045
confint(Pro.CGREC.pt)
## Computing profile confidence intervals ...
## 2.5 % 97.5 %
## .sig01 0.3918434 0.9030728
## .sig02 0.0000000 0.6318090
## .sigma 1.8002131 1.9942528
## (Intercept) 11.6337283 12.7277087
## PatchU -3.4089867 -2.6718668
## TreatmentSpring Only -0.6397023 0.6672886
fixef(Pro.CGREC.pt)
## (Intercept) PatchU TreatmentSpring Only
## 12.188124924 -3.040368577 0.006071581
summary(Pro.CGREC.p)
## Linear mixed model fit by maximum likelihood ['lmerMod']
## Formula: PROTEIN ~ Patch + (1 | Pasture/Month)
## Data: subset(ProGO, REC == "CGREC")
##
## AIC BIC logLik deviance df.resid
## 3160.0 3183.1 -1575.0 3150.0 752
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.1774 -0.6825 -0.1494 0.5193 4.9335
##
## Random effects:
## Groups Name Variance Std.Dev.
## Month:Pasture (Intercept) 0.377738 0.61460
## Pasture (Intercept) 0.008784 0.09372
## Residual 3.583837 1.89310
## Number of obs: 757, groups: Month:Pasture, 24; Pasture, 8
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 12.1913 0.2148 56.75
## PatchU -3.0405 0.1875 -16.22
##
## Correlation of Fixed Effects:
## (Intr)
## PatchU -0.729
confint(Pro.CGREC.p) #p or pt give basically the same CI for protein; pt accounts for # treatment
## Computing profile confidence intervals ...
## 2.5 % 97.5 %
## .sig01 0.3918435 0.9030232
## .sig02 0.0000000 0.6317193
## .sigma 1.8002166 1.9942563
## (Intercept) 11.7528694 12.6227237
## PatchU -3.4087572 -2.6724256
#Pro.CGREC.F <- lmer(PROTEIN ~ Patch + Treatment + (1|Month) + (1|Pasture),
# data=subset(ProGO, REC=="CGREC"), REML=FALSE)
# summary(Pro.CGREC.F)
# confint(Pro.CGREC.F)
Biomass
Bio.HREC.null <- lmer(Biomass ~ 1 + (1|Pasture/Month),
data=subset(ProGO, REC=="HREC"), REML=FALSE)
Bio.HREC.p <- lmer(Biomass ~ Patch + (1|Pasture/Month),
data=subset(ProGO, REC=="HREC"), REML=FALSE)
Bio.HREC.pt <- lmer(Biomass ~ Patch + Treatment + (1|Pasture/Month),
data=subset(ProGO, REC=="HREC"), REML=FALSE)
anova(Bio.HREC.null, Bio.HREC.p, Bio.HREC.pt) %>%
pander(caption="comparison of the three models")
comparison of the three models (continued below)
| Bio.HREC.null |
4 |
3011 |
3027 |
-1501 |
3003 |
NA |
NA |
| Bio.HREC.p |
5 |
2974 |
2994 |
-1482 |
2964 |
39.05 |
1 |
| Bio.HREC.pt |
6 |
2976 |
2999 |
-1482 |
2964 |
0.2868 |
1 |
| Bio.HREC.null |
NA |
| Bio.HREC.p |
4.139e-10 |
| Bio.HREC.pt |
0.5923 |
#.P AND PT ARE ABOUT THE SAME
anova(Bio.HREC.null, Bio.HREC.p) %>%
pander(caption="null vs patch")
null vs patch (continued below)
| Bio.HREC.null |
4 |
3011 |
3027 |
-1501 |
3003 |
NA |
NA |
| Bio.HREC.p |
5 |
2974 |
2994 |
-1482 |
2964 |
39.05 |
1 |
| Bio.HREC.null |
NA |
| Bio.HREC.p |
4.139e-10 |
anova(Bio.HREC.null, Bio.HREC.p, Bio.HREC.pt) %>%
pander(caption="null vs patch and treatment")
null vs patch and treatment (continued below)
| Bio.HREC.null |
4 |
3011 |
3027 |
-1501 |
3003 |
NA |
NA |
| Bio.HREC.p |
5 |
2974 |
2994 |
-1482 |
2964 |
39.05 |
1 |
| Bio.HREC.pt |
6 |
2976 |
2999 |
-1482 |
2964 |
0.2868 |
1 |
| Bio.HREC.null |
NA |
| Bio.HREC.p |
4.139e-10 |
| Bio.HREC.pt |
0.5923 |
Bio.HREC.INT <- lmer(Biomass ~ Patch * Treatment + (1|Pasture/Month),
data=subset(ProGO, REC=="HREC"), REML=FALSE)
anova(Bio.HREC.pt, Bio.HREC.INT)
## Data: subset(ProGO, REC == "HREC")
## Models:
## Bio.HREC.pt: Biomass ~ Patch + Treatment + (1 | Pasture/Month)
## Bio.HREC.INT: Biomass ~ Patch * Treatment + (1 | Pasture/Month)
## Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq)
## Bio.HREC.pt 6 2975.6 2999.3 -1481.8 2963.6
## Bio.HREC.INT 7 2967.0 2994.7 -1476.5 2953.0 10.565 1 0.001152 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary(Bio.HREC.pt)
## Linear mixed model fit by maximum likelihood ['lmerMod']
## Formula: Biomass ~ Patch + Treatment + (1 | Pasture/Month)
## Data: subset(ProGO, REC == "HREC")
##
## AIC BIC logLik deviance df.resid
## 2975.6 2999.3 -1481.8 2963.6 381
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.8403 -0.7404 -0.1928 0.5383 2.9856
##
## Random effects:
## Groups Name Variance Std.Dev.
## Month:Pasture (Intercept) 4.832 2.198
## Pasture (Intercept) 15.985 3.998
## Residual 117.168 10.824
## Number of obs: 387, groups: Month:Pasture, 18; Pasture, 6
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 15.572 2.721 5.723
## PatchU 8.250 1.286 6.416
## TreatmentSheep 1.949 3.598 0.542
##
## Correlation of Fixed Effects:
## (Intr) PatchU
## PatchU -0.357
## TreatmntShp -0.660 0.000
confint(Bio.HREC.pt) #PATCH SIG, TREATMENT NOT
## Computing profile confidence intervals ...
## 2.5 % 97.5 %
## .sig01 0.000000 4.522105
## .sig02 1.760247 8.627888
## .sigma 10.088122 11.654854
## (Intercept) 9.450430 21.690614
## PatchU 5.723598 10.777187
## TreatmentSheep -6.375216 10.303160
summary(Bio.HREC.p)
## Linear mixed model fit by maximum likelihood ['lmerMod']
## Formula: Biomass ~ Patch + (1 | Pasture/Month)
## Data: subset(ProGO, REC == "HREC")
##
## AIC BIC logLik deviance df.resid
## 2973.9 2993.7 -1481.9 2963.9 382
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.8305 -0.7365 -0.1879 0.5464 2.9937
##
## Random effects:
## Groups Name Variance Std.Dev.
## Month:Pasture (Intercept) 4.832 2.198
## Pasture (Intercept) 16.898 4.111
## Residual 117.172 10.825
## Number of obs: 387, groups: Month:Pasture, 18; Pasture, 6
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 16.545 2.082 7.948
## PatchU 8.250 1.286 6.416
##
## Correlation of Fixed Effects:
## (Intr)
## PatchU -0.467
confint(Bio.HREC.p) #PATCH IS SIGNIFICANT
## Computing profile confidence intervals ...
## 2.5 % 97.5 %
## .sig01 0.000000 4.522122
## .sig02 1.867843 8.843722
## .sigma 10.088254 11.655059
## (Intercept) 11.969645 21.132414
## PatchU 5.723341 10.777057
Bio.CGREC.null <- lmer(Biomass ~ 1 + (1|Pasture/Month),
data=subset(ProGO, REC=="CGREC"), REML=FALSE)
Bio.CGREC.p <- lmer(Biomass ~ Patch + (1|Pasture/Month),
data=subset(ProGO, REC=="CGREC"), REML=FALSE)
Bio.CGREC.pt <- lmer(Biomass ~ Patch + Treatment + (1|Pasture/Month),
data=subset(ProGO, REC=="CGREC"), REML=FALSE)
anova(Bio.CGREC.null, Bio.CGREC.p, Bio.CGREC.pt) %>%
pander(caption="comparison of the three models")
comparison of the three models (continued below)
| Bio.CGREC.null |
4 |
5499 |
5517 |
-2745 |
5491 |
NA |
NA |
| Bio.CGREC.p |
5 |
5429 |
5453 |
-2710 |
5419 |
71.18 |
1 |
| Bio.CGREC.pt |
6 |
5430 |
5458 |
-2709 |
5418 |
1.298 |
1 |
| Bio.CGREC.null |
NA |
| Bio.CGREC.p |
3.255e-17 |
| Bio.CGREC.pt |
0.2546 |
#.P IS THE BETTER MODEL
anova(Bio.CGREC.null, Bio.CGREC.p) %>%
pander(caption="null vs patch")
null vs patch (continued below)
| Bio.CGREC.null |
4 |
5499 |
5517 |
-2745 |
5491 |
NA |
NA |
| Bio.CGREC.p |
5 |
5429 |
5453 |
-2710 |
5419 |
71.18 |
1 |
| Bio.CGREC.null |
NA |
| Bio.CGREC.p |
3.255e-17 |
anova(Bio.CGREC.null, Bio.CGREC.pt) %>%
pander(caption="null vs patch and treatment")
null vs patch and treatment (continued below)
| Bio.CGREC.null |
4 |
5499 |
5517 |
-2745 |
5491 |
NA |
NA |
| Bio.CGREC.pt |
6 |
5430 |
5458 |
-2709 |
5418 |
72.48 |
2 |
| Bio.CGREC.null |
NA |
| Bio.CGREC.pt |
1.823e-16 |
Bio.CGREC.INT <- lmer(Biomass ~ Patch * Treatment + (1|Pasture/Month),
data=subset(ProGO, REC=="CGREC"), REML=FALSE)
anova(Bio.CGREC.pt, Bio.CGREC.INT)
## Data: subset(ProGO, REC == "CGREC")
## Models:
## Bio.CGREC.pt: Biomass ~ Patch + Treatment + (1 | Pasture/Month)
## Bio.CGREC.INT: Biomass ~ Patch * Treatment + (1 | Pasture/Month)
## Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq)
## Bio.CGREC.pt 6 5430.2 5457.9 -2709.1 5418.2
## Bio.CGREC.INT 7 5427.3 5459.7 -2706.7 5413.3 4.8385 1 0.02783 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary(Bio.CGREC.pt)
## Linear mixed model fit by maximum likelihood ['lmerMod']
## Formula: Biomass ~ Patch + Treatment + (1 | Pasture/Month)
## Data: subset(ProGO, REC == "CGREC")
##
## AIC BIC logLik deviance df.resid
## 5430.2 5457.9 -2709.1 5418.2 751
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.1457 -0.6682 -0.1832 0.4933 4.8393
##
## Random effects:
## Groups Name Variance Std.Dev.
## Month:Pasture (Intercept) 5.292e+00 2.301e+00
## Pasture (Intercept) 2.978e-13 5.457e-07
## Residual 7.237e+01 8.507e+00
## Number of obs: 757, groups: Month:Pasture, 24; Pasture, 8
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 15.4692 1.0799 14.325
## PatchU 7.3321 0.8429 8.698
## TreatmentSpring Only 1.3027 1.1283 1.155
##
## Correlation of Fixed Effects:
## (Intr) PatchU
## PatchU -0.679
## TrtmntSprnO -0.552 0.053
confint(Bio.CGREC.pt)
## Computing profile confidence intervals ...
## Warning in nextpar(mat, cc, i, delta, lowcut, upcut): Last two rows have
## identical or NA .zeta values: using minstep
## Warning in nextpar(mat, cc, i, delta, lowcut, upcut): Last two rows have
## identical or NA .zeta values: using minstep
## Warning in FUN(X[[i]], ...): non-monotonic profile for .sig02
## Warning in confint.thpr(pp, level = level, zeta = zeta): bad spline fit
## for .sig02: falling back to linear interpolation
## 2.5 % 97.5 %
## .sig01 1.438861 3.448834
## .sig02 0.000000 2.144707
## .sigma 8.089442 8.961712
## (Intercept) 13.316869 17.611006
## PatchU 5.677675 8.986246
## TreatmentSpring Only -1.057132 3.662826
fixef(Bio.CGREC.pt)
## (Intercept) PatchU TreatmentSpring Only
## 15.469173 7.332050 1.302719
summary(Bio.CGREC.p)
## Linear mixed model fit by maximum likelihood ['lmerMod']
## Formula: Biomass ~ Patch + (1 | Pasture/Month)
## Data: subset(ProGO, REC == "CGREC")
##
## AIC BIC logLik deviance df.resid
## 5429.5 5452.6 -2709.7 5419.5 752
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.1674 -0.6667 -0.1902 0.5015 4.8265
##
## Random effects:
## Groups Name Variance Std.Dev.
## Month:Pasture (Intercept) 5.712 2.390
## Pasture (Intercept) 0.000 0.000
## Residual 72.370 8.507
## Number of obs: 757, groups: Month:Pasture, 24; Pasture, 8
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 16.1534 0.9105 17.741
## PatchU 7.2843 0.8419 8.652
##
## Correlation of Fixed Effects:
## (Intr)
## PatchU -0.772
confint(Bio.CGREC.p)
## Computing profile confidence intervals ...
## 2.5 % 97.5 %
## .sig01 1.455832 3.561241
## .sig02 0.000000 2.395668
## .sigma 8.089484 8.961780
## (Intercept) 14.325825 17.954553
## PatchU 5.631343 8.936701
Veg Structure
CI for patch difference
kable(VOR, caption="Overall structure difference and 95% CI")
Overall structure difference and 95% CI
| CGREC |
Mixed-Native |
Spring Only |
-33 |
-14 |
-22.9 |
| CGREC |
Mixed-Native |
Spring + Summer |
-49 |
-10 |
-29.5 |
| HREC |
Mixed-Introduced |
Cattle |
-74 |
-12 |
-43.4 |
| HREC |
Mixed-Introduced |
Sheep |
-46 |
-39 |
-42.5 |
ggplot(VOR, aes(x=Treatment, Y=sim.mean, color=Treatment))+
geom_point(aes(y=sim.mean, shape=Location),size=5)+
scale_color_brewer(palette="Set1")+
scale_fill_brewer(palette="Set1")+
geom_hline(yintercept = 0, color="black", linetype="dashed")+
geom_errorbar(aes(ymin=ci.L, ymax=ci.U),size=2, width=0.2)+
coord_cartesian(ylim=c(-100,20))+
scale_y_continuous(breaks=c(-80,-60,-40,-20,0,20))+
xlab ("Treatment")+
ylab ("Vegetative Structure (% Difference) + 95 % CI")+
guides(shape=guide_legend(title="Location:"), color=guide_legend(title="Treatment:"))+
theme_bw()+
theme (axis.text=element_text(size=12, face="bold"),
#legend.position = c(.5, .88),
legend.position="bottom",
panel.grid.minor = element_blank(),
panel.grid.major = element_blank(),
axis.title=element_text(size=14, face="bold"),
plot.title=element_text(size=12, face="bold", hjust=0.5))
