I am trying to analyse some data and was given R code to do this with but
there seem to be errors in the code. My level of knowledge is improving but
still limited.
The details are;
Data on clover lines; Lines.txt attached. Comma seperations
Code:
options(digits=3)
clover <- read.table("Lines.txt",header=T,sep=",")
vnames <- names(clover);nv <- length(vnames)
flevels <- levels(clover$Line)
par(oma=c(0,6,4,6),mfrow=c(2,1),mar=c(4,4,1,0))
clover$Y <- clover[,"LeafLmm"]
mn.tab <- tapply(clover$Y,list(clover$Line),mean)
sd.tab <- sqrt(tapply(clover$Y,list(clover$Line),var))
om <- order(mn.tab)
plot(mn.tab,sd.tab,las=1,xlab="mean",ylab="s.d.")
lines(lowess(mn.tab,sd.tab))
model.i <- glm(Y ~ C(clover$Line,Cline),data=clover,family=gaussian),
Beta <- summary(model.i)$coefficients
new.order <- order(Beta[-1,1]) + 1
Obeta <- Beta[c(1,new.order),]
detect.diff <- quantile(Obeta,probs=seq(0.2,0.8,0.2))
std.dev <- (min(Obeta[-1,2]) + max(Obeta[-1,2]))/2*sqrt(5)
sink("ClovAttr.txt")
print(paste("_____________",vnames[(i+2)],"_____________"))
cat("__________ mean table _______ \n")
print(mn.tab[om])
cat("__________ s.d. table _______ \n")
print(sd.tab[om])
print(model.i$family)
print(Obeta)
print(detect.diff)
cat("_________________________________\n")
sink()
library(ctest)
stdev <- round(std.dev*c(1.25,1,0.75),2)
stderr <- round(stdev/sqrt(5),2)
samp.size <- 3:10
beta <- 0.8
ylab <- "detectable differences"
for (j in seq(along=stdev)){
Delta <- numeric(0)
for(k in seq(along=samp.size)){
Power.calcs <- power.t.test(n=samp.size[k],delta=NULL,sd=stdev[j],
sig.level=0.05,power=beta,type="one.sample",alternative="one.sided")
beta <- Power.calcs$power
Delta <- c(Delta,Power.calcs$delta)
} # end of the k loop
if(j==1){
max.Delta <- max(Delta)
plot(samp.size,Delta,type='l',xlab="n",ylab=ylab,las=1,
ylim=c(0,max.Delta),xlim=c(3,10))
} # end of if(j==1)
else lines(samp.size,Delta,lty=j)
} # end of j loop
mtext(vnames[i],side=3,outer=T,cex=1.2)
legend(7,max.Delta,legend=paste("s.e",stderr),lty=1:length(stdev))
text(5,0.9*max.Delta,"power=0.8")
_________________________________
The issue is the number of leaves to be sampled.
As above I only needed to analyse as a linear model, Plot curves of
detectable difference @ 3 standard deviations etc. All the code is there
but there seems to be errors. model.i I don't know what Cline is supposed
to be (possibly clover$Line), the detect.diff line is incorrect, the looping
through J&K appears incorrect and I don't know what i in the vnames is
supposed to be.
Any help? This is probably simple stuff for experienced people.
Stephen Beale
-------------- next part --------------
An embedded and charset-unspecified text was scrubbed...
Name: Lines.txt
Url:
https://stat.ethz.ch/pipermail/r-help/attachments/20030521/20fd52eb/Lines.txt