Post date: Dec 17, 2013 8:13:56 PM
I am trying an alternative, and probably better approach to estimate breeding values for wild SLA and GLA individuals. Basically, I want to fit a single BSLMM probit model for survival on Ms or Ac with experimental and wild SLA and GLA individuals. I want to fit the model directly to the binary response. I think this is important as "correcting" for differences in survival between source populations will likely "correct" differences in estimated breeding values in the wild (this so-called correction removes the potential actual effect). Of course the downside of this new approach is that there will be excess LD (mainly captured through the polygenic term via the kinship matrix). This will perhaps allow us to do a better job of explaining variation in survival than we would otherwise, but I think this is ok in the current context. Besides, applying this model to all populations is really the more important test.
First I need to format the data.
I used estpost to generate genotype mean files (ghat_* in lycaeides_hostplant/genotypes/).
I ran the R script formatGenotypes.R to create genotype matrixes (gmat_*) with for gemma. There is one matrix per host plant treatment and each includes genotypes for experimental and wild individuals. This script also generates survival phenotype files with individuals in the same order as the genotype data and with NA's for wild-caught individuals.
I then ran makeGemmaIn.pl (still in genotypes folder) to add the locus information to the genotype data and make geno_* files for gemma.
I moved the geno and pheno files to melGemma and ran gemma for posterior prediction with the probit model, 5 million iterations, 1 million burnin, thinning interval of 40, and 5 chains for each host plant. These are in the long queue on the dorc.
These jobs were just sitting in the long queue. I re-started them in the batch queue (96 hours) with slightly different MCMC conditions, 10 chains, 3 million iterations, 1 million iteration burnin, and thinning interval of 40. THe new job ids are 68299-68313.