BasicASReml: Difference between revisions

From RAJ INFO
Jump to navigation Jump to search
Created page with "Category:R <source lang=R> # set working directory to "course material" setwd("H:\CORP\ITCRD-CORPRND\Agri_Common\GenStat_ASReml_VSNi\ICT DAY2\ICT DAY2 PART1\MET") lib..."
 
mNo edit summary
Line 83: Line 83:
diag.asr <- update(diag.asr)
diag.asr <- update(diag.asr)
summary(diag.asr)$varcomp
summary(diag.asr)$varcomp
----
# set working directory to "course material"
setwd("R:/CORP/ITCRD-CORPRND/Agri_Common/GenStat_ASReml_VSNi/ICT DAY2/ICT DAY2 PART1/ITC")
getwd()
rm(list = ls())
ALLD = read.csv("ALL.csv", header=TRUE)
summary(ALLD)
hist(ALLD$gBH48)
is.factor(ALLD$Genotype)
is.factor(ALLD$LOC)
#ALLD$GENO = as.factor(ALLD$Genotype)
#ALLD = ALLD[,-14]
library(asreml)
# fixed effect for Genotype + LOC + Genotype.LOC
ALLD.asr1 <- asreml(gBH48 ~ Genotype*LOC, data=ALLD)
str(ALLD.asr1)
ALLD.asr1$vcoeff
ALLD.asr1$coefficients
effects.fix = coef(ALLD.asr1)$fixed
# loc as random, genotype as fixed #
# fixed effect for Genotype + LOC + Genotype.LOC
ALLD.asr3 <- asreml(gBH48 ~ Genotype, random = ~ LOC + Genotype.LOC, data=ALLD)
#str(ALLD.asr3)
effects.fix3 = coef(ALLD.asr3)$fixed
write.csv(effects.fix3,"effects.csv")
# random effect for Genotype + LOC + Genotype.LOC
ALLD.asr2 <- asreml(gBH48 ~ 1, random = ~ Genotype+LOC+Genotype:LOC, data=ALLD)
summary(ALLD.asr2)$varcomp
#h2 = 4.5/(7.26+4.5+14.58)
bv<-coef(ALLD.asr2)$random
ALLD.asr4 <- asreml(gBH48 ~ Genotype, random = ~ LOC + Genotype:LOC, data=ALLD)
#str(ALLD.asr4)
effects.fix4 = coef(ALLD.asr4)$fixed
write.csv(effects.fix4,"effects4.csv")


</source>
</source>

Revision as of 09:20, 19 November 2016


# set working directory to "course material"
setwd("H:\CORP\ITCRD-CORPRND\Agri_Common\GenStat_ASReml_VSNi\ICT DAY2\ICT DAY2 PART1\MET")


library(asreml)
#library(myf)

# read data
sbt <- asreml.read.table("sbt_all.dat",sep="",header=TRUE,na.strings="*")

# sort data frame
sbt <- sbt[order(sbt$Trial,sbt$Row,sbt$Col),]

# get cycle factor for T3
sbt$colcycle <- as.factor(((as.numeric(sbt$Col)-1)%%4))
# get cycle factor for T1
sbt$rowcycle <- as.factor(((as.numeric(sbt$Row)-1)%%2))

# calculate variate of interest
sbt$dmyield <- sbt$yield*sbt$dmpc/100

# get subsets for individual trials
sbt.T1 <- subset(sbt,Trial=="T1")
sbt.T2 <- subset(sbt,Trial=="T2")
sbt.T3 <- subset(sbt,Trial=="T3")
sbt.T4 <- subset(sbt,Trial=="T4")
sbt.T5 <- subset(sbt,Trial=="T5")

# establish models for individual trials 

# T1
T1.asr <- asreml(fixed=dmyield~1+rowcycle, random=~Variety+Rep+Col+units,
                  rcov=~ar1(Row):ar1(Col), na.method.X="include", data=sbt.T1)
summary(T1.asr)$varcomp
metplot(T1.asr)
plot(T1.asr)
anova(T1.asr)

# T2
T2.asr <- asreml(fixed=dmyield~1, random=~Variety+Rep+Col+units,
                 rcov=~ar1(Row):ar1(Col), na.method.X="include", data=sbt.T2)
T2.asr <- update(T2.asr)
summary(T2.asr)$varcomp
plot(T2.asr)
metplot(T2.asr)

# T3
T3.asr <- asreml(fixed=dmyield~1+colcycle, random=~Variety+Rep+Col+units,
                 rcov=~ar1(Row):ar1(Col), na.method.X="include", data=sbt.T3)
summary(T3.asr)$varcomp
plot(T3.asr)
metplot(T3.asr)

# T4
T4.asr <- asreml(fixed=dmyield~1, random=~Variety+Rep,
                 rcov=~ar1(Row):ar1(Col), na.method.X="include", data=sbt.T4)
summary(T4.asr)$varcomp
plot(T4.asr)
metplot(T4.asr)

# T5
T5.asr <- asreml(fixed=dmyield~1+lin(Row), random=~Variety+Rep,
                 rcov=~id(Row):ar1(Col), na.method.X="include", data=sbt.T5)
summary(T5.asr)$varcomp
plot(T5.asr)
#metplot(T5.asr)

# one-stage joint analysis: independence across trials for GxE

diag.asr <- asreml(fixed=dmyield~1+Trial+at(Trial,1):rowcycle+
                         at(Trial,3):colcycle+at(Trial,5):lin(Row),
                   random=~diag(Trial):Variety + at(Trial):Rep
                          + at(Trial,c(1:3)):Col+ at(Trial,c(1:3)):units,
                   rcov=~at(Trial,c(1:4)):ar1(Row):ar1(Col)
                        +at(Trial,5):id(Row):ar1(Col),
                   na.method.X="include", data=sbt)
# use update function to get convergence
diag.asr <- update(diag.asr)
summary(diag.asr)$varcomp

----
# set working directory to "course material"
setwd("R:/CORP/ITCRD-CORPRND/Agri_Common/GenStat_ASReml_VSNi/ICT DAY2/ICT DAY2 PART1/ITC")
getwd()

rm(list = ls())

ALLD = read.csv("ALL.csv", header=TRUE)
summary(ALLD)

hist(ALLD$gBH48)
is.factor(ALLD$Genotype)
is.factor(ALLD$LOC)

#ALLD$GENO = as.factor(ALLD$Genotype)
#ALLD = ALLD[,-14]

library(asreml)

# fixed effect for Genotype + LOC + Genotype.LOC 
ALLD.asr1 <- asreml(gBH48 ~ Genotype*LOC, data=ALLD)
str(ALLD.asr1)
ALLD.asr1$vcoeff
ALLD.asr1$coefficients

effects.fix = coef(ALLD.asr1)$fixed

# loc as random, genotype as fixed #
# fixed effect for Genotype + LOC + Genotype.LOC 
ALLD.asr3 <- asreml(gBH48 ~ Genotype, random = ~ LOC + Genotype.LOC, data=ALLD)
#str(ALLD.asr3)

effects.fix3 = coef(ALLD.asr3)$fixed
write.csv(effects.fix3,"effects.csv")

# random effect for Genotype + LOC + Genotype.LOC 
ALLD.asr2 <- asreml(gBH48 ~ 1, random = ~ Genotype+LOC+Genotype:LOC, data=ALLD)
summary(ALLD.asr2)$varcomp
#h2 = 4.5/(7.26+4.5+14.58)
bv<-coef(ALLD.asr2)$random


ALLD.asr4 <- asreml(gBH48 ~ Genotype, random = ~ LOC + Genotype:LOC, data=ALLD)
#str(ALLD.asr4)
effects.fix4 = coef(ALLD.asr4)$fixed
write.csv(effects.fix4,"effects4.csv")