CTFSworkshop(data=TRUE)
?rnbinom
?rnbinom
rnbinom(1,size=10,mu=.2)
rnbinom(1,size=10,mu=.2)
?rnbinom
rnbinom(1,size=1,mu=.2)
rnbinom(25,size=1,mu=.2)
rnbinom(25,size=.1,mu=.2)
rnbinom(25,size=100,mu=.2)
rnbinom(25,size=100,mu=.8)
r=rnbinom(25,size=100,mu=.8)
mean(r)
r=rnbinom(2500,size=100,mu=.8)
mean(r)
table(r)
rnorm(10,mean=15,sd=3)
source('~/meetings_workshops/Rworkshops/teaching.functions.r')
source('~/meetings_workshops/Rmodeling/teaching.functions.r')
sample.2level
sample.2level(populations=5,popN=10,overallmean=5,overallSD=1,withinSD=5)
source('~/meetings_workshops/Rmodeling/teaching.functions.r')
sample.2level(populations=5,popN=10,overallmean=5,overallSD=1,withinSD=5)
sample.2level(populations=5,popN=10,overallmean=5,overallSD=0,withinSD=5)
source('~/meetings_workshops/Rmodeling/teaching.functions.r')
sample.2level(populations=5,popN=10,overallmean=5,overallSD=0,withinSD=5)
popmean
overallsd
overallSD
n
y
popN
withinN
withinSD
popmean[i]
observedmean
n
observedmean
n
c
Q
source('~/meetings_workshops/Rmodeling/teaching.functions.r')
sample.2level(populations=5,popN=10,overallmean=5,overallSD=0,withinSD=5)
result=sample.2level(populations=100,popN=10,overallmean=5,overallSD=0,withinSD=5)
head(result)
source('~/meetings_workshops/Rmodeling/teaching.functions.r')
source('~/meetings_workshops/Rmodeling/teaching.functions.r')
debug(sample.2level)
result=sample.2level(populations=100,popN=10,overallmean=5,overallSD=0,withinSD=5)
c
Q
undebug(sample.2level)
result=sample.2level(populations=100,popN=10,overallmean=5,overallSD=0,withinSD=5)
head(result)
hist(result$observedmean)
sd(result$observedmean)
result=sample.2level(populations=100,popN=10,overallmean=5,overallSD=0,withinSD=15)
head(result)
graphics.off()
par(mfcol=c(2,1))
result=sample.2level(populations=100,popN=10,overallmean=5,overallSD=0,withinSD=5)
hist(result$observedmean)
result=sample.2level(populations=100,popN=10,overallmean=5,overallSD=0,withinSD=15)
hist(result$observedmean)
sd(result$observedmean)
result=sample.2level(populations=100,popN=10,overallmean=5,overallSD=0,withinSD=5)
sd(result$observedmean)
hist(result$observedmean)
result=sample.2level(populations=100,popN=10,overallmean=5,overallSD=2,withinSD=5)
hist(result$observedmean)
sd(result$observedmean)
result=sample.2level(populations=100,popN=1000,overallmean=5,overallSD=2,withinSD=5)
sd(result$observedmean)
hist(result$observedmean)
?lmer
search()
options(with=60)
options(width=60)
search()
ls(6)
head(grwfull300)
options(width=80)
head(grwfull300)
options(width=70)
head(grwfull300)
ls(6)
search()
head(grwfull300)
colnames(grwfull300)
library(lme4)
?lmer
mean(subset(grwfull300,Species=='beilpe')$incgr)
sd(subset(grwfull300,Species=='beilpe')$incgr)
spmean=tapply(grwfull300$incgr,grwfull300$Species,mean)
spSD=tapply(grwfull300$incgr,grwfull300$Species,sd)
length(spmean)
head(spmean)
head(spmean,10)
head(spSD,10)
mean(grwfull300$incgr)
mean(spmean)
hist(spmean)
graphics.off()
hist(spmean)
which(spmean>15)
which(spmean>20)
spmean[which(spmean>20)]
dim(subset(grwfull300,sp=='tremmi'))
dim(subset(grwfull300,Species=='tremmi'))
dim(subset(grwfull300,Species=='pit1ma'))
spmean=lm(incgr~Species,data=grwfull300)
head(coef(spmean))
spmean=tapply(grwfull300$incgr,grwfull300$Species,mean)
hist(subset(grwfull300,Species=='alsebl')$incgr)
mod=lmer(incgr~1+(1|Species),data=grwfull300)
fit=lm(y~0+x)
mod=lmer(incgr~1+(1|Species),data=grwfull300)
fixef(mod)
mean(spmean)
mean(grwfull300$incgr)
head(coef(mod)$Species)
str(mod$coef)
head(spmean)
dim(subset(grwfull300,Species=='amaico'))
(subset(grwfull300,Species=='amaico'))
plot(spmean,coef(mod)$Species)
plot(spmean,coef(mod)$Species[,1])
head(coef(mod)$Species)
abline(0,1)
plot(spmean,coef(mod)$Species[,1])
identify(spmean,coef(mod)$Species[,1],names(spmean))
head(spmean)
spmean=tapply(grwfull300$incgr,grwfull300$Species,mean)
head(names(spmean))
head(names(spmean))
identify(spmean,coef(mod)$Species[,1],names(spmean))
(subset(grwfull300,Species=='guetfo'))
coef(mod)$Species['guetfo',]
head(coef(mod)$Species)
abline(0,1)
dim(subset(grwfull300,Species=='tachve'))
head(spmean)
display(mod)
raresp=names(which(abund<10))
abund=table(grwfull300$Species)
raresp=names(which(abund<10))
abundsp=names(which(abund>50))
plot(as.numeric(abund),spmean)
plot(spmean,coef(mod)$Species[,1])
identify(spmean,coef(mod)$Species[,1],names(spmean))
head(abund)
plot(as.numeric(abund),spmean)
length(raresp)
abundsp
raresp
sd(spmean)
sd(spmean[abundsp])
sd(spmean[raresp])
mean(spmean)
mean(spmean[abundsp])
growth=seq(.1,20,by=.1)
prob=dnorm(growth,mean=mean(spmean[abundsp]),sd=sd(spmean[abundsp]))
plot(growth,prob,type='l')
plot(spmean,coef(mod)$Species[,1])
plot(as.numeric(abund),spmean)
plot(growth,prob,type='l')
mod=lmer(incgr~1+(1|Species),data=grwfull300)
display(mod)
fixef(mod)
display(mod)
colnames(grwfull300)
modLn=lmer(LnGrowth~1+(1|Species),data=grwfull300)
display(modLn)
head(coef(modLn)$Species)
sd(coef(modLn)$Species[,1])
mean(coef(modLn)$Species[,1])
hist(coef(modLn)$Species[,1])
hist(coef(modLn)$Species[,1])
hist(grwfull300$incgr)
hist(grwfull300$LnGrowth)
modLn=lmer(LnGrowth~1+(1|Species),data=grwfull300)
head(se.coef(modLn)$Species)
head(coef(modLn)$Species)
modLn
spmean=tapply(grwfull300$LnGrowth,grwfull300$Species,mean)
head(spmean)
head(coef(modLn)$Species)
head(abund)
colnames(grwfull300)
lm(LnGrowth~dbh1,data=grwfull300)
summary(lm(LnGrowth~dbh1,data=grwfull300))
modLn=lmer(LnGrowth~1+(1|Species),data=grwfull300)
lm(LnGrowth~1+dbh1,data=grwfull300)
lm(LnGrowth~1+dbh1+Species+dbh1*Species,data=grwfull300)
modReg1=lmer(LnGrowth~1+dbh1+(1|Species),data=grwfull300)
head(coef(modReg1)$Species)
fixef(modReg1)
modReg2=lmer(LnGrowth~1+dbh1+(1+dbh1|Species),data=grwfull300)
modReg2=lmer(LnGrowth~1+dbh1+(1+dbh1|Species),data=grwfull300)
head(coef(modReg2)$Species)
head(abund)
fixef(modReg1)
fixef(modReg2)
modReg1
display(modReg1)
display(modReg1,digits=5)
summary(lm(LnGrowth~dbh1,data=grwfull300))
display(modReg1,digits=5)
display(modReg2,digits=5)
head(coef(modReg2)$Species)
(coef(modReg2)$Species$dbh1)
hist(coef(modReg2)$Species$dbh1)
hist(coef(modReg2)$Species$dbh1,breaks=50)
hist(coef(modReg2)$Species$dbh1,breaks=25)
plot(coef(modReg2)$Species)
colnames(grwfull300)
grwfull300$dbhctr=grwfull300$dbh1-mean(grwfull300$dbh1)
plot(LnGrowth~dbh1,data=grwfull300)
plot(LnGrowth~dbhctr,data=grwfull300)
modReg3=lmer(LnGrowth~1+dbhctr+(1+dbhctr|Species),data=grwfull300)
display(modReg3)
display(modReg3,digits=5)
display(modReg2,digits=5)
display(modReg3,digits=5)
modReg4=lmer(LnGrowth~1+dbhctr+(1|Species),data=grwfull300)
head(coef(modReg4)$Species)
modReg5=lmer(LnGrowth~1+dbhctr+(0+dbhctr|Species),data=grwfull300)
head(coef(modReg5)$Species)
AIC(modReg2,modReg3,modReg4,modReg5)
anova(modReg2,modReg3,modReg4,modReg5)
cf=as.matrix(coef(modReg3)$Species)
dim(cf)
head(cf)
plot(LnGrowth~dbh1,data=grwfull300)
plot(LnGrowth~dbhctr,data=grwfull300)
plot(LnGrowth~dbhctr,data=grwfull300,cex=.5,pch=16)
abline(cf[1,])
for(i in 1:146) abline(cf[i,],col='gray')
points(LnGrowth~dbhctr,data=grwfull300,cex=.5,pch=16)
for(i in 1:146) abline(cf[i,],col='orange')
fixef(modReg3)
abline(fixef(modReg3),col='black',lwd=2)
subset(cf,cf[,1]>2)
points(LnGrowth~dbhctr,data=subset(grwfull300,Species=='tachve'),cex=.5,pch=16)
points(LnGrowth~dbhctr,data=subset(grwfull300,Species=='tachve'),cex=.5,pch=16,col='red')
fit=lm(LnGrowth~1,data=grwfull300)
summary(fit)
anova(fit)
aic(fit)
AIC(fit)
AIC(modReg5,fit)
AIC(modReg5,modReg3,fit)
savehistory('~/meetings_workshops/Rmodeling/Cenpat/history11Oct.txt')