In fboehm/pleio: Likelihood ratio test to distinguish pleiotropy from linkage

The linear mixed effects model

$$ vec(Y) = XB + vec(G) + vec(E)$$ where $Y$ is a n by 2 matrix of phenotypes, with one row per subject and one column per trait, $X$ is a 2n by 2f block diagonal matrix consisting of two n by f matrices, one for each of two loci. $n$ is the number of subjects and $f$ is the number of founders. $G$ is a n by 2 matrix of random effects, with $$G \sim MN(0, K, V_g)$$ where $MN(0, K, V_g)$ denotes the matrix-variate normal distribution with mean zero, row covariance $K$ and column covariance $V_g$. $E$ is a n by 2 matrix of random errors with one row per subject and one trait per column. We assume that $$E \sim MN(0, I_n, V_e)$$

Load DO Data

library(pleiotropy)
library(dplyr)
# load packages
library(broman) # contains the winsorize() function
library(tidyverse)
# load data
### genotype probabilities ("probs") in form used by R/qtl2
### see ../R/0_DOQTLprobs2qtl2.R for how they were converted
PATH_TO_DERIVED_DATA <- "~/Box Sync/attie/attiedo"
PATH_TO_DERIVED_DATA <- "~/attie"
load(file.path(PATH_TO_DERIVED_DATA, "DerivedData/GM_Attie_allele_call_haploprobs_4qtl2_wave5.Rdata"))
### clinical phenotypes + phenotype dictionary
### ("pheno_clin" and "pheno_clin_dict")
load(file.path(PATH_TO_DERIVED_DATA, "DerivedData/pheno_clin.RData"))
### kinship matrices ("loco" method) ("K")
load(file.path(PATH_TO_DERIVED_DATA, "DerivedData/kinship.RData"))
### covariate matrix ("covar")
#load(file.path(PATH_TO_DERIVED_DATA, "DerivedData/covar.RData"))
### physical map of the markers in the probs array
load(file.path(PATH_TO_DERIVED_DATA, "DerivedData/probs_pmap.RData"))