Description Format Details Source References Examples

Wheat varieties with heteroskedastic yields

A data frame with 52 observations on the following 3 variables.

`env`

environment, 13 levels

`gen`

genotype, 4 levels

`yield`

yield

Yield of 4 varieties of wheat at 13 locations in Oklahoma, USA.

The data was used to explore variability between varieties.

F. A. Graybill, 1954.
Variance heterogeneity in a randomized block design,
*Biometrics*, 10, 516-520.

Hans-Pieter Piepho, 1994.
Missing observations in the analysis of stability.
*Heredity*, 72, 141–145.
http://doi.org/10.1038/hdy.1994.20

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 | ```
data(graybill.heteroskedastic)
dat <- graybill.heteroskedastic
# Genotypes are obviously not homoscedastic
boxplot(yield ~ gen, dat, main="graybill.heteroskedastic")
# Shukla stability variance of each genotype, same as Grubbs' estimate
# Matches Piepho 1994 page 143.
# Do not do this! Nowadays, use mixed models instead.
if(require("reshape2")){
datm <- acast(dat, gen~env)
w <- datm
w <- sweep(w, 1, rowMeans(datm))
w <- sweep(w, 2, colMeans(datm))
w <- w + mean(datm)
w <- rowSums(w^2)
k=4; n=13
sig2 <- k*w/((k-2)*(n-1)) - sum(w)/((k-1)*(k-2)*(n-1))
## sig2
## G1 G2 G3 G4
## 145.98 -14.14 75.15 18.25
var.shukla <- function(x,N){
# Estimate variance of shukla stability statistics
# Piepho 1994 equation (5)
K <- length(x) # num genotypes
S <- outer(x,x)
S1 <- diag(S)
S2 <- rowSums(S) - S1
S[!upper.tri(S)] <- 0 # Make S upper triangular
# The ith element of S3 is the sum of the upper triangular elements of S,
# excluding the ith row and ith column
S3 <- sum(S) - rowSums(S) - colSums(S)
var.si2 <- 2*S1/(N-1) + 4/( (N-1)*(K-1)^2 ) * ( S2 + S3/(K-2)^2 )
return(var.si2)
}
# Set negative estimates to zero
sig2[sig2<0] <- 0
# Variance of shukla stat. Match Piepho 1994, table 5, example 1
var.shukla(sig2,13)
## G1 G2 G3 G4
## 4069.3296 138.9424 1423.0797 306.5270
}
``` |

agridat documentation built on Nov. 30, 2017, 1:02 a.m.

Embedding an R snippet on your website

Add the following code to your website.

For more information on customizing the embed code, read Embedding Snippets.