%\VignetteIndexEntry{Application Tutorial: OmicKriging} 
\documentclass[a4paper]{article}

\title{Application Tutorial: OmicKriging}
\author{Keston Aquino-Michaels, Heather E. Wheeler, Vassily V. Trubetskoy and Hae Kyung Im}


\usepackage{hyperref}
\usepackage{Sweave} 

\SweaveOpts{engine=R} 

\begin{document}
\SweaveOpts{concordance=TRUE}
\maketitle
Method citation: Wheeler HE, Aquino-Michaels K, Gamazon ER, Trubetskoy VV, Dolan ME, Huang RS, Cox NJ, Im HK. Poly-Omic Prediction of Complex Traits: OmicKriging. Genet Epidemiol. 2014 May 2. doi: 10.1002/gepi.21808. 
\begin{center}
  \line(1,0){250}
\end{center}
\begin{center}
  \line(1,0){250}
\end{center}

\begin{section}{Running OmicKriging with Example Data}

To install from CRAN:
\begin{Schunk}
\begin{Sinput}
> install.packages("OmicKriging")
\end{Sinput}
\end{Schunk}
\begin{center}
\line(1,0){250}
\end{center}
Start by loading OmicKriging functions into R:
\begin{Schunk}
\begin{Sinput}
> library(OmicKriging)
\end{Sinput}
\end{Schunk}
Define paths to the genetic relatedness (GCTA binary format), gene expression, and phenotype data files (paths may differ based on where the files are located). The path.package() function returns the package installation directory. These files will later be passed to upcoming functions:
\begin{Schunk}
<<>>=
library(OmicKriging)

binaryFile <- system.file(package = "OmicKriging",
                "doc/vignette_data/ig_genotypes.grm.bin")
binaryFileBase <- substr(binaryFile,1, nchar(binaryFile) - 4)
expressionFile <- system.file(package = "OmicKriging",
                    "doc/vignette_data/ig_gene_subset.txt.gz") 
phenotypeFile <- system.file(package = "OmicKriging",
                  "doc/vignette_data/ig_pheno.txt") 
@

\end{Schunk}
Load the phenotype data into R:
\begin{Schunk}
<<>>=
pheno <- read.table(phenotypeFile, header = T)
@

\end{Schunk}
Load a pre-computed GCTA GRM into R:
\begin{Schunk}
<<>>=
grmMat <- read_GRMBin(binaryFileBase)
@

\end{Schunk}
By default, grmFilePrefix is set to NULL. However, if specified, this function will save the computed GRM to disk in GCTA binary format. Additionally by default both snpList and sampleList are set to NULL. However you may restrict the GRM calculation by specifying a vector of sample IDs or a vector of SNP IDs.
\\
\\
Load and calculate a gene expression relatedness matrix (GXM) with the following function:
\begin{Schunk}
<<>>=
gxmMat <- make_GXM(expFile = expressionFile)
@

\end{Schunk}
Similarly, by default, gxmFilePrefix is set to NULL, however if specified, this function will save the computed GXM to disk in GCTA binary format.
\begin{center}
\line(1,0){250}
\end{center}
Additional convenience functions are included to perform principal components analysis (PCA):
\begin{Schunk}
<<>>=
pcMatXM <- make_PCs_irlba(gxmMat, n.top = 10)

pcMatGM <- make_PCs_irlba(grmMat, n.top = 10)

pcMat <- cbind(pcMatGM, pcMatXM[match(rownames(pcMatGM), rownames(pcMatXM)),])
@

\end{Schunk}
\begin{center}
\line(1,0){250}
\end{center}
The following convenience function allows the user to perform n-fold cross-validation. Specify the number of cores you wish to use (default = "all"), the number of cross-validation folds desired (default = 10), covariates (by default covar.mat = NULL), the phenotype object, pheno.id (by default = 1 (the first phenotype in the file)), the h2 vector and a list of the correlation matrices to be included.
 \\
\\
Note: The sum of the h2 vector must be between 0 and 1. In this example, we will give each matrix equal weight.

\begin{Schunk}

<<>>=
result <- krigr_cross_validation(pheno.df = pheno,
	cor.list = list(grmMat, gxmMat),
	h2.vec = c(0.5, 0.5),
	covar.mat = pcMat,
	ncore = 2,
	nfold = "LOOCV")
@

\end{Schunk}
This function will return a data.frame with column Ypred corresponding to the predicted values and column Ytest corresponding to the measured phenotypes.
\begin{center}
\line(1,0){250}
\end{center}
\begin{center}
Congratulations! 
\\
You have just completed the OmicKriging tutorial!
\end{center}
\begin{center}
\line(1,0){250}
\end{center}

\end{section}

\begin{section}{Cleanup}

This is a cleanup step for the vignette on Windows; typically not
needed for users.

\begin{Schunk}
<<closeConnetions,results=hide>>=
allCon <- showConnections()
socketCon <- as.integer(rownames(allCon)[allCon[, "class"] == "sockconn"])
sapply(socketCon, function(ii) close.connection(getConnection(ii)) )
@
\end{Schunk}
\end{section}

\end{document}