Calculate Heritability
H2.fun <- function(x){
pdata <- x
pdata = split(pdata$pheno, pdata$isotype)
pdata.notNAcnt = sapply(pdata, function(x){sum(!is.na(x))})
pdata[pdata.notNAcnt<2]=NULL
pdata.melted = melt(pdata)
names(pdata.melted)=c('pheno', 'isotype')
pdata.melted$strain=as.factor(pdata.melted$isotype)
reffMod = lme4::lmer(pheno ~ 1 + (1|strain), data=pdata.melted)
Var_Random_effect <- as.numeric(VarCorr(reffMod))
Var_Residual <- attr(VarCorr(reffMod), "sc")^2
H2 <- Var_Random_effect/(Var_Random_effect+Var_Residual)
H2
}