HW4
Hongyu Chen
Monday, October 13, 2014
This is an R Markdown document. Markdown is a simple formatting syntax for authoring HTML, PDF, and MS Word documents. For more details on using R Markdown see http://rmarkdown.rstudio.com.
When you click the Knit button a document will be generated that includes both content as well as the output of any embedded R code chunks within the document. You can embed an R code chunk like this:
Fertilization at planting influences seedling growth
Hongyu Chen
RPI, RIN:661405156
Oct.11 V1.0
1. Setting
System under test
Data of this test is obtained from an article about influence of fertilization on plants. The aim of this analysis is to test data given in this paper and try to reach a same conclusion as the author’s. Data is from Table2 provided in the article, about response of foliar nutrient concentrations to fertilization treatments. A quick view of this dataset is as below.
data<-read.csv("fertilizer.csv", header=TRUE)
data## Species Fertilizer Nitrogen Phosphate Potassium
## 1 Aspen SR-H 3.30 0.15 1.10
## 2 Aspen SR-L 2.98 0.15 0.99
## 3 Aspen FR-H 3.22 0.16 1.08
## 4 Aspen FR-L 3.12 0.15 1.12
## 5 Aspen IA-H 2.87 0.16 1.10
## 6 Aspen IA-L 2.36 0.14 0.88
## 7 Aspen UC 2.21 0.17 0.79
## 8 White Spruce SR-H 1.72 0.08 0.54
## 9 White Spruce SR-L 1.52 0.08 0.55
## 10 White Spruce FR-H 1.72 0.08 0.61
## 11 White Spruce FR-L 1.51 0.08 0.53
## 12 White Spruce IA-H 1.37 0.08 0.60
## 13 White Spruce IA-L 1.03 0.10 0.56
## 14 White Spruce UC 0.85 0.08 0.60summary(data)## Species Fertilizer Nitrogen Phosphate
## Aspen :7 FR-H:2 Min. :0.85 Min. :0.080
## White Spruce:7 FR-L:2 1st Qu.:1.51 1st Qu.:0.080
## IA-H:2 Median :1.97 Median :0.120
## IA-L:2 Mean :2.13 Mean :0.119
## SR-H:2 3rd Qu.:2.95 3rd Qu.:0.150
## SR-L:2 Max. :3.30 Max. :0.170
## UC :2
## Potassium
## Min. :0.530
## 1st Qu.:0.570
## Median :0.700
## Mean :0.789
## 3rd Qu.:1.058
## Max. :1.120
## Factors and Levels
There are two factors: species and fertilization treatment. In the factor of ‘Species’ there are two species as different levels; in factor of ‘Fertilizer’, there are seven different types of fertilizers as different levels. This is a two factor multiple level model.
Continuous variables (if any)
Nutrition concentrations are continuous variables, in specific are concentration of nitrogen, phosphate and potassium.
Response variables
Concentration of different nutrients are set as response variables, which are ‘Nitrogen’, ‘Phosphate’ and ‘Potassium’ in the data set.
The Data: How is it organized and what does it look like?
Data used in this test is obtained from a published paper, Table 2. Original table is about influence of fertilizer treatment on concentration of nutrients in different species seedlings. Two factors are species and fertilizer treatments.
Randomization
For each experiment a randomized complete block design was employed.
2. (Experimental) Design
How will the experiment be organized and conducted to test the hypothesis?
Read data from csv document, first go through exploratory analysis, then use ANOVA to test the effect of each factor and the combination on response variables. Use Tukey’s multiple pairwise comparison to identify significant differences between fertilizer treatments and model checking methods to check adequacy of a selected model as example. Through that we I hope to examine effect from both factors respectively.
What is the rationale for this design?
Fertilization may alleviate plants nutrient deficiencies but broadcast fertilization with immediately available fertilizers (IAF) results in generally low rates of nutrient recovery for planted trees. Directed application of controlled-release fertilizer (CRF) to the rhizosphere offers an alternative to extend nutrient longevity while reducing nutrient leaching or uptake by competing vegetation. Test about effect of fertilizers is conducted with white spruce and aspen.
Randomize: What is the Randomization Scheme?
For each experiment, a randomized complete block design was employed. Specifically, plant density, numbers of plant population in each block and possible influence from plants nearby are blocked.
Replicate: Are there replicates, blocking and/or repeated measures?
There were four replicate blocks in this design. Each block contained 84 experimental seedlings hand-planted at a spacing of 2m X 2m with each block surrounded by a buffer row of non-experimental seedlings planted at the same density. Also repeated measurement was applied.
3. (Statistical) Analysis
(Exploratory Data Analysis) Graphics and descriptive summary
# Logical vector identifying all factors
spe<-factor(data$Species)
fert<-factor(data$Fertilizer)
Nitrogen<-factor(data$Nitrogen)
# Quick view of samples collected
levels(spe)## [1] "Aspen" "White Spruce"levels(fert)## [1] "FR-H" "FR-L" "IA-H" "IA-L" "SR-H" "SR-L" "UC"boxplot(Nitrogen~spe, data=data)
boxplot(Nitrogen~fert, data=data)
boxplot(Phosphate~spe, data=data)
boxplot(Phosphate~fert, data=data)
boxplot(Potassium~spe, data=data)
boxplot(Potassium~fert, data=data)
From boxplots above, different species show an obviously different nutrients concentration, therefore species is probably a factor to explain variance in nutrient concentration, including nitrogen, phosphate and potassium. However, with different fertilizer treatments, nitrogen and potassium concentration vary relatively great, but difference in terms of phosphate median is not that obvious, indicating that fertilizer may be able to partly explain variance of nutrient concentration in plants seedlings.
Testing
ANOVA
#variance of species
modela1=aov(Nitrogen ~ spe, data=data)
anova(modela1)## Analysis of Variance Table
##
## Response: Nitrogen
## Df Sum Sq Mean Sq F value Pr(>F)
## spe 1 7.64 7.64 52.4 1e-05 ***
## Residuals 12 1.75 0.15
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1modela2=aov(Phosphate ~ spe, data=data)
anova(modela2)## Analysis of Variance Table
##
## Response: Phosphate
## Df Sum Sq Mean Sq F value Pr(>F)
## spe 1 0.01786 0.01786 234 3.1e-09 ***
## Residuals 12 0.00091 0.00008
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1modela3=aov(Potassium ~ spe, data=data)
anova(modela3)## Analysis of Variance Table
##
## Response: Potassium
## Df Sum Sq Mean Sq F value Pr(>F)
## spe 1 0.673 0.673 76.7 1.5e-06 ***
## Residuals 12 0.105 0.009
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1All P-values are much smaller than 0.05, which means there is very small probability that variance in nutrient concentration with regards to species is due to ramdomization.
#variance of fertilizer
modelb1=aov(Nitrogen ~ fert, data=data)
anova(modelb1)## Analysis of Variance Table
##
## Response: Nitrogen
## Df Sum Sq Mean Sq F value Pr(>F)
## fert 6 1.72 0.286 0.26 0.94
## Residuals 7 7.67 1.096modelb2=aov(Phosphate ~ fert, data=data)
anova(modelb2)## Analysis of Variance Table
##
## Response: Phosphate
## Df Sum Sq Mean Sq F value Pr(>F)
## fert 6 0.00017 0.000029 0.01 1
## Residuals 7 0.01860 0.002657modelb3=aov(Potassium ~ fert, data=data)
anova(modelb3)## Analysis of Variance Table
##
## Response: Potassium
## Df Sum Sq Mean Sq F value Pr(>F)
## fert 6 0.046 0.0077 0.07 1
## Residuals 7 0.732 0.1046P-values are greater than 0.05, which means basically under such experiment settings, nutrient concentration is independent to fertilizer treatments. Probably each type of fertilizer is able to provide similar amount of nutrients for plants, which may play an important role in plant growth.
#variance of interaction of both factors
modelc1=aov(Nitrogen ~ spe*fert, data=data)
anova(modelc1)## Warning: ANOVA F-tests on an essentially perfect fit are unreliable## Analysis of Variance Table
##
## Response: Nitrogen
## Df Sum Sq Mean Sq F value Pr(>F)
## spe 1 7.64 7.64
## fert 6 1.72 0.29
## spe:fert 6 0.03 0.01
## Residuals 0 0.00modelc2=aov(Phosphate ~ spe*fert, data=data)
anova(modelc2)## Warning: ANOVA F-tests on an essentially perfect fit are unreliable## Analysis of Variance Table
##
## Response: Phosphate
## Df Sum Sq Mean Sq F value Pr(>F)
## spe 1 0.01786 0.01786
## fert 6 0.00017 0.00003
## spe:fert 6 0.00074 0.00012
## Residuals 0 0.00000modelc3=aov(Potassium ~ spe*fert, data=data)
anova(modelc3)## Warning: ANOVA F-tests on an essentially perfect fit are unreliable## Analysis of Variance Table
##
## Response: Potassium
## Df Sum Sq Mean Sq F value Pr(>F)
## spe 1 0.673 0.673
## fert 6 0.046 0.008
## spe:fert 6 0.059 0.010
## Residuals 0 0.000It can be seen that for each type of nutrient, interaction of species and fertilizer treatments leads to a very small P-vale, also there is a warning message: ANOVA F-tests on an essentially perfect fit are unreliable, indicating interaction of species and fertilizer is a nearly perfect model to explain variance of nutrient concentration.
Tukey’s Honest Significant Difference Test
TukeyHSD(modelc1, ordered = FALSE, conf.level = 0.95)## Warning: 产生了NaNs
## Warning: 产生了NaNs
## Warning: 产生了NaNs## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = Nitrogen ~ spe * fert, data = data)
##
## $spe
## diff lwr upr p adj
## White Spruce-Aspen -1.477 NaN NaN NaN
##
## $fert
## diff lwr upr p adj
## FR-L-FR-H -0.155 NaN NaN NaN
## IA-H-FR-H -0.350 NaN NaN NaN
## IA-L-FR-H -0.775 NaN NaN NaN
## SR-H-FR-H 0.040 NaN NaN NaN
## SR-L-FR-H -0.220 NaN NaN NaN
## UC-FR-H -0.940 NaN NaN NaN
## IA-H-FR-L -0.195 NaN NaN NaN
## IA-L-FR-L -0.620 NaN NaN NaN
## SR-H-FR-L 0.195 NaN NaN NaN
## SR-L-FR-L -0.065 NaN NaN NaN
## UC-FR-L -0.785 NaN NaN NaN
## IA-L-IA-H -0.425 NaN NaN NaN
## SR-H-IA-H 0.390 NaN NaN NaN
## SR-L-IA-H 0.130 NaN NaN NaN
## UC-IA-H -0.590 NaN NaN NaN
## SR-H-IA-L 0.815 NaN NaN NaN
## SR-L-IA-L 0.555 NaN NaN NaN
## UC-IA-L -0.165 NaN NaN NaN
## SR-L-SR-H -0.260 NaN NaN NaN
## UC-SR-H -0.980 NaN NaN NaN
## UC-SR-L -0.720 NaN NaN NaN
##
## $`spe:fert`
## diff lwr upr p adj
## White Spruce:FR-H-Aspen:FR-H -1.500e+00 NaN NaN NaN
## Aspen:FR-L-Aspen:FR-H -1.000e-01 NaN NaN NaN
## White Spruce:FR-L-Aspen:FR-H -1.710e+00 NaN NaN NaN
## Aspen:IA-H-Aspen:FR-H -3.500e-01 NaN NaN NaN
## White Spruce:IA-H-Aspen:FR-H -1.850e+00 NaN NaN NaN
## Aspen:IA-L-Aspen:FR-H -8.600e-01 NaN NaN NaN
## White Spruce:IA-L-Aspen:FR-H -2.190e+00 NaN NaN NaN
## Aspen:SR-H-Aspen:FR-H 8.000e-02 NaN NaN NaN
## White Spruce:SR-H-Aspen:FR-H -1.500e+00 NaN NaN NaN
## Aspen:SR-L-Aspen:FR-H -2.400e-01 NaN NaN NaN
## White Spruce:SR-L-Aspen:FR-H -1.700e+00 NaN NaN NaN
## Aspen:UC-Aspen:FR-H -1.010e+00 NaN NaN NaN
## White Spruce:UC-Aspen:FR-H -2.370e+00 NaN NaN NaN
## Aspen:FR-L-White Spruce:FR-H 1.400e+00 NaN NaN NaN
## White Spruce:FR-L-White Spruce:FR-H -2.100e-01 NaN NaN NaN
## Aspen:IA-H-White Spruce:FR-H 1.150e+00 NaN NaN NaN
## White Spruce:IA-H-White Spruce:FR-H -3.500e-01 NaN NaN NaN
## Aspen:IA-L-White Spruce:FR-H 6.400e-01 NaN NaN NaN
## White Spruce:IA-L-White Spruce:FR-H -6.900e-01 NaN NaN NaN
## Aspen:SR-H-White Spruce:FR-H 1.580e+00 NaN NaN NaN
## White Spruce:SR-H-White Spruce:FR-H -4.441e-16 NaN NaN NaN
## Aspen:SR-L-White Spruce:FR-H 1.260e+00 NaN NaN NaN
## White Spruce:SR-L-White Spruce:FR-H -2.000e-01 NaN NaN NaN
## Aspen:UC-White Spruce:FR-H 4.900e-01 NaN NaN NaN
## White Spruce:UC-White Spruce:FR-H -8.700e-01 NaN NaN NaN
## White Spruce:FR-L-Aspen:FR-L -1.610e+00 NaN NaN NaN
## Aspen:IA-H-Aspen:FR-L -2.500e-01 NaN NaN NaN
## White Spruce:IA-H-Aspen:FR-L -1.750e+00 NaN NaN NaN
## Aspen:IA-L-Aspen:FR-L -7.600e-01 NaN NaN NaN
## White Spruce:IA-L-Aspen:FR-L -2.090e+00 NaN NaN NaN
## Aspen:SR-H-Aspen:FR-L 1.800e-01 NaN NaN NaN
## White Spruce:SR-H-Aspen:FR-L -1.400e+00 NaN NaN NaN
## Aspen:SR-L-Aspen:FR-L -1.400e-01 NaN NaN NaN
## White Spruce:SR-L-Aspen:FR-L -1.600e+00 NaN NaN NaN
## Aspen:UC-Aspen:FR-L -9.100e-01 NaN NaN NaN
## White Spruce:UC-Aspen:FR-L -2.270e+00 NaN NaN NaN
## Aspen:IA-H-White Spruce:FR-L 1.360e+00 NaN NaN NaN
## White Spruce:IA-H-White Spruce:FR-L -1.400e-01 NaN NaN NaN
## Aspen:IA-L-White Spruce:FR-L 8.500e-01 NaN NaN NaN
## White Spruce:IA-L-White Spruce:FR-L -4.800e-01 NaN NaN NaN
## Aspen:SR-H-White Spruce:FR-L 1.790e+00 NaN NaN NaN
## White Spruce:SR-H-White Spruce:FR-L 2.100e-01 NaN NaN NaN
## Aspen:SR-L-White Spruce:FR-L 1.470e+00 NaN NaN NaN
## White Spruce:SR-L-White Spruce:FR-L 1.000e-02 NaN NaN NaN
## Aspen:UC-White Spruce:FR-L 7.000e-01 NaN NaN NaN
## White Spruce:UC-White Spruce:FR-L -6.600e-01 NaN NaN NaN
## White Spruce:IA-H-Aspen:IA-H -1.500e+00 NaN NaN NaN
## Aspen:IA-L-Aspen:IA-H -5.100e-01 NaN NaN NaN
## White Spruce:IA-L-Aspen:IA-H -1.840e+00 NaN NaN NaN
## Aspen:SR-H-Aspen:IA-H 4.300e-01 NaN NaN NaN
## White Spruce:SR-H-Aspen:IA-H -1.150e+00 NaN NaN NaN
## Aspen:SR-L-Aspen:IA-H 1.100e-01 NaN NaN NaN
## White Spruce:SR-L-Aspen:IA-H -1.350e+00 NaN NaN NaN
## Aspen:UC-Aspen:IA-H -6.600e-01 NaN NaN NaN
## White Spruce:UC-Aspen:IA-H -2.020e+00 NaN NaN NaN
## Aspen:IA-L-White Spruce:IA-H 9.900e-01 NaN NaN NaN
## White Spruce:IA-L-White Spruce:IA-H -3.400e-01 NaN NaN NaN
## Aspen:SR-H-White Spruce:IA-H 1.930e+00 NaN NaN NaN
## White Spruce:SR-H-White Spruce:IA-H 3.500e-01 NaN NaN NaN
## Aspen:SR-L-White Spruce:IA-H 1.610e+00 NaN NaN NaN
## White Spruce:SR-L-White Spruce:IA-H 1.500e-01 NaN NaN NaN
## Aspen:UC-White Spruce:IA-H 8.400e-01 NaN NaN NaN
## White Spruce:UC-White Spruce:IA-H -5.200e-01 NaN NaN NaN
## White Spruce:IA-L-Aspen:IA-L -1.330e+00 NaN NaN NaN
## Aspen:SR-H-Aspen:IA-L 9.400e-01 NaN NaN NaN
## White Spruce:SR-H-Aspen:IA-L -6.400e-01 NaN NaN NaN
## Aspen:SR-L-Aspen:IA-L 6.200e-01 NaN NaN NaN
## White Spruce:SR-L-Aspen:IA-L -8.400e-01 NaN NaN NaN
## Aspen:UC-Aspen:IA-L -1.500e-01 NaN NaN NaN
## White Spruce:UC-Aspen:IA-L -1.510e+00 NaN NaN NaN
## Aspen:SR-H-White Spruce:IA-L 2.270e+00 NaN NaN NaN
## White Spruce:SR-H-White Spruce:IA-L 6.900e-01 NaN NaN NaN
## Aspen:SR-L-White Spruce:IA-L 1.950e+00 NaN NaN NaN
## White Spruce:SR-L-White Spruce:IA-L 4.900e-01 NaN NaN NaN
## Aspen:UC-White Spruce:IA-L 1.180e+00 NaN NaN NaN
## White Spruce:UC-White Spruce:IA-L -1.800e-01 NaN NaN NaN
## White Spruce:SR-H-Aspen:SR-H -1.580e+00 NaN NaN NaN
## Aspen:SR-L-Aspen:SR-H -3.200e-01 NaN NaN NaN
## White Spruce:SR-L-Aspen:SR-H -1.780e+00 NaN NaN NaN
## Aspen:UC-Aspen:SR-H -1.090e+00 NaN NaN NaN
## White Spruce:UC-Aspen:SR-H -2.450e+00 NaN NaN NaN
## Aspen:SR-L-White Spruce:SR-H 1.260e+00 NaN NaN NaN
## White Spruce:SR-L-White Spruce:SR-H -2.000e-01 NaN NaN NaN
## Aspen:UC-White Spruce:SR-H 4.900e-01 NaN NaN NaN
## White Spruce:UC-White Spruce:SR-H -8.700e-01 NaN NaN NaN
## White Spruce:SR-L-Aspen:SR-L -1.460e+00 NaN NaN NaN
## Aspen:UC-Aspen:SR-L -7.700e-01 NaN NaN NaN
## White Spruce:UC-Aspen:SR-L -2.130e+00 NaN NaN NaN
## Aspen:UC-White Spruce:SR-L 6.900e-01 NaN NaN NaN
## White Spruce:UC-White Spruce:SR-L -6.700e-01 NaN NaN NaN
## White Spruce:UC-Aspen:UC -1.360e+00 NaN NaN NaNDataset check
data## Species Fertilizer Nitrogen Phosphate Potassium
## 1 Aspen SR-H 3.30 0.15 1.10
## 2 Aspen SR-L 2.98 0.15 0.99
## 3 Aspen FR-H 3.22 0.16 1.08
## 4 Aspen FR-L 3.12 0.15 1.12
## 5 Aspen IA-H 2.87 0.16 1.10
## 6 Aspen IA-L 2.36 0.14 0.88
## 7 Aspen UC 2.21 0.17 0.79
## 8 White Spruce SR-H 1.72 0.08 0.54
## 9 White Spruce SR-L 1.52 0.08 0.55
## 10 White Spruce FR-H 1.72 0.08 0.61
## 11 White Spruce FR-L 1.51 0.08 0.53
## 12 White Spruce IA-H 1.37 0.08 0.60
## 13 White Spruce IA-L 1.03 0.10 0.56
## 14 White Spruce UC 0.85 0.08 0.60Dataset is complete without any missing information, perpaps Tukey’s HSD test cannot be applied to essentially perfit fitted ANOVA. Then try to analyze species and fertilizer respectively, take nitrogen as an example.
Tukey_species<-TukeyHSD(modela1, ordered = FALSE, conf.level = 0.95)
Tukey_species## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = Nitrogen ~ spe, data = data)
##
## $spe
## diff lwr upr p adj
## White Spruce-Aspen -1.477 -1.922 -1.033 0Tukey_fertilizer<-TukeyHSD(modelb1, ordered = FALSE, conf.level = 0.95)
Tukey_fertilizer## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = Nitrogen ~ fert, data = data)
##
## $fert
## diff lwr upr p adj
## FR-L-FR-H -0.155 -4.304 3.994 1.0000
## IA-H-FR-H -0.350 -4.499 3.799 0.9998
## IA-L-FR-H -0.775 -4.924 3.374 0.9844
## SR-H-FR-H 0.040 -4.109 4.189 1.0000
## SR-L-FR-H -0.220 -4.369 3.929 1.0000
## UC-FR-H -0.940 -5.089 3.209 0.9615
## IA-H-FR-L -0.195 -4.344 3.954 1.0000
## IA-L-FR-L -0.620 -4.769 3.529 0.9950
## SR-H-FR-L 0.195 -3.954 4.344 1.0000
## SR-L-FR-L -0.065 -4.214 4.084 1.0000
## UC-FR-L -0.785 -4.934 3.364 0.9834
## IA-L-IA-H -0.425 -4.574 3.724 0.9994
## SR-H-IA-H 0.390 -3.759 4.539 0.9996
## SR-L-IA-H 0.130 -4.019 4.279 1.0000
## UC-IA-H -0.590 -4.739 3.559 0.9961
## SR-H-IA-L 0.815 -3.334 4.964 0.9802
## SR-L-IA-L 0.555 -3.594 4.704 0.9972
## UC-IA-L -0.165 -4.314 3.984 1.0000
## SR-L-SR-H -0.260 -4.409 3.889 1.0000
## UC-SR-H -0.980 -5.129 3.169 0.9538
## UC-SR-L -0.720 -4.869 3.429 0.9892plot(Tukey_species)
plot(Tukey_fertilizer)
Tukey’s HSD test can tell us which groups in the sample differ significantly. For individual test if there is a p-adj.value<0.05, then there is a statistical difference between the mean response variables of those two levels, which is because of something other than randomization. Therefore for Tukey’s test of species, p-adj value is smaller than 0.05, demonstrating there is difference between different species levels caused by something other than randomization. However for Tukey’s test of fertilizer, all p-adj values are larger than 0.05, indicating there is no statistical difference between fertilizer samples, and variance is attributed to randomization. This is a similar conclusion as from ANOVA.
Diagnostics/Model Adequacy Checking
Take nitrogen concentration as an example. Set the interaction of species and fertilizer treatments as the model.
modelc1=aov(Nitrogen~spe*fert, data=data)
qqnorm(residuals(modelc1))
qqline(residuals(modelc1))
qqnorm plot and qqline of residuals exhibit perfect linear pattern of residuals which are all falling at 0, indicating an essentially perfect model.
interaction.plot(data$Species,data$Fertilizer,data$Nitrogen)
Interaction polt does not appear as parallel lines, which means the interaction model is probably adequate to explain variance of nitrogen concentration.
plot(fitted(modelc1),residuals(modelc1))
Plot of fitted model and residuals model shows residuals are exactly distributed at 0, indicating a perfect model.
4. References to the literature
5. Appendices
A summary of, or pointer to, the raw data
Fertilization at planting influences seedling growth and vegetative competition on a post-mining boreal reclamation site Joshua L. Sloan, Douglass F. Jacobs http://link.springer.com/article/10.1007/s11056-013-9378-4