Nothing
projectR Vignette
In projectR: Functions for the projection of weights from PCA, CoGAPS, NMF, correlation, and clustering
knitr::opts_chunk$set(echo = TRUE)
options(scipen = 1, digits = 2)
set.seed(1234)
options(tinytex.verbose = TRUE)
Introduction
Technological advances continue to spur the exponential growth of biological data as illustrated by the rise of the omics—genomics, transcriptomics, epigenomics, proteomics, etc.—each with there own high throughput technologies. In order to leverage the full power of these resources, methods to integrate multiple data sets and data types must be developed. The reciprocal nature of the genomic, transcriptomic, epigenomic, and proteomic biology requires that the data provides a complementary view of cellular function and regulatory organization; however, the technical heterogeneity and massive size of high-throughput data even within a particular omic makes integrated analysis challenging. To address these challenges, we developed projectR, an R package for integrated analysis of high dimensional omic data. projectR uses the relationships defined within a given high dimensional data set, to interrogate related biological phenomena in an entirely new data set. By relying on relative comparisons within data type, projectR is able to circumvent many issues arising from technological variation. For a more extensive example of how the tools in the projectR package can be used for in silico experiments, or additional information on the algorithm, see Stein-O'Brien, et al.
Getting started with projectR
Installation Instructions
For automatic Bioconductor package installation, start R, and run:
BiocManager::install("projectR")
Methods
Projection can roughly be defined as a mapping or transformation of points from one space to another often lower dimensional space. Mathematically, this can described as a function $\varphi(x)=y : \Re^{D} \mapsto \Re^{d}$ s.t. $d \leq D$ for $x \in \Re^{D}, y \in \Re^{d}$ @Barbakh:2009bw . The projectR package uses projection functions defined in a training dataset to interrogate related biological phenomena in an entirely new data set. These functions can be the product of any one of several methods common to "omic" analyses including regression, PCA, NMF, clustering. Individual sections focussing on one specific method are included in the vignette. However, the general design of the projectR function is the same regardless.
The base projectR function
The generic projectR function is executed as follows:
library(projectR)
projectR(data, loadings, dataNames=NULL, loadingsNames=NULL, NP = NULL, full = false)
Input Arguments
The inputs that must be set each time are only the data and loadings, with all other inputs having default values. However, incongruities in the feature mapping between the data and loadings, i.e. a different format for the rownames of each object, will throw errors or result in an empty mapping and should be checked before running. To overcoming mismatched feature names in the objects themselve, the /code{dataNames} and /code{loadingNames} arguements can be manually supplied by the user.
The arguments are as follows:
\begin{description}
\item[data]{a dataset to be projected into the pattern space}
\item[loadings]{a matrix of continous values with unique rownames to be projected}
\item[dataNames]{a vector containing unique name, i.e. gene names, for the rows of the target dataset to be used to match features with the loadings, if not provided by \texttt{rownames(data)}. Order of names in vector must match order of rows in data.}
\item[loadingsNames]{a vector containing unique names, i.e. gene names, for the rows of loadings to be used to match features with the data, if not provided by \texttt{rownames(loadings)}. Order of names in vector must match order of rows in loadings.}
\item[NP]{vector of integers indicating which columns of loadings object to use. The default of NP = NA will use entire matrix.}
\item[full]{logical indicating whether to return the full model solution. By default only the new pattern object is returned.}
\end{description}
The \texttt{loadings} arguement in the generic projectR function is suitable for use with any genernal feature space, or set of feature spaces, whose rows annotation links them to the data to be projected. Ex: the coeffients associated with individual genes as the result of regression analysis or the amplituded values of individual genes as the result of non-negative matrix factorization (NMF).
Output
The basic output of the base projectR function, i.e. \texttt{full=FALSE}, returns \texttt{projectionPatterns} representing relative weights for the samples from the new data in this previously defined feature space, or set of feature spaces. The full output of the base projectR function, i.e. \texttt{full=TRUE}, returns \texttt{projectionFit}, a list containing \texttt{projectionPatterns} and \texttt{Projection}. The \texttt{Projection} object contains additional information from the proceedure used to obtain the \texttt{projectionPatterns}. For the the the base projectR function, \texttt{Projection} is the full lmFit model from the package r BiocStyle::Biocpkg("limma")
PCA projection
Projection of principle components is achieved by matrix multiplication of a new data set by previously generated eigenvectors, or gene loadings. If the original data were standardized such that each gene is centered to zero average expression level, the principal components are normalized eigenvectors of the covariance matrix of the genes. Each PC is ordered according to how much of the variation present in the data they contain. Projection of the original samples into each PC will maximized the variance of the samples in the direction of that component and uncorrelated to previous components. Projection of new data places the new samples into the PCs defined by the original data. Because the components define an orthonormal basis set, they provide an isomorphism between a vector space, $V$, and $\Re^n$ which preserves inner products. If $V$ is an inner product space over $\Re$ with orthonormal basis $B = v_1,...,v_n$ and $v \epsilon V s.t [v]_B = (r_1,...,r_n)$, then finding the coordinate of $v_i$ in $v$ is precisely the inner product of $v$ with $v_i$, i.e. $r_i = \langle v,v_i \rangle$. This formulation is implemented for only those genes belonging to both the new data and the PC space. The \texttt{projectR} function has S4 method for class \texttt{prcomp}.
Obtaining PCs to project.
# data to define PCs
library(projectR)
data(p.RNAseq6l3c3t)
# do PCA on RNAseq6l3c3t expression data
pc.RNAseq6l3c3t<-prcomp(t(p.RNAseq6l3c3t))
pcVAR <- round(((pc.RNAseq6l3c3t$sdev)^2/sum(pc.RNAseq6l3c3t$sdev^2))*100,2)
dPCA <- data.frame(cbind(pc.RNAseq6l3c3t$x,pd.RNAseq6l3c3t))
#plot pca
library(ggplot2)
setCOL <- scale_colour_manual(values = c("blue","black","red"), name="Condition:")
setFILL <- scale_fill_manual(values = c("blue","black","red"),guide = FALSE)
setPCH <- scale_shape_manual(values=c(23,22,25,25,21,24),name="Cell Line:")
pPCA <- ggplot(dPCA, aes(x=PC1, y=PC2, colour=ID.cond, shape=ID.line,
fill=ID.cond)) +
geom_point(aes(size=days),alpha=.6)+
setCOL + setPCH + setFILL +
scale_size_area(breaks = c(2,4,6), name="Day") +
theme(legend.position=c(0,0), legend.justification=c(0,0),
legend.direction = "horizontal",
panel.background = element_rect(fill = "white",colour=NA),
legend.background = element_rect(fill = "transparent",colour=NA),
plot.title = element_text(vjust = 0,hjust=0,face="bold")) +
labs(title = "PCA of hPSC PolyA RNAseq",
x=paste("PC1 (",pcVAR[1],"% of varience)",sep=""),
y=paste("PC2 (",pcVAR[2],"% of varience)",sep=""))
Projecting prcomp objects
# data to project into PCs from RNAseq6 l3c3t expression data
data(p.ESepiGen4c1l)
library(projectR)
PCA2ESepi <- projectR(data = p.ESepiGen4c1l$mRNA.Seq,loadings=pc.RNAseq6l3c3t,
full=TRUE, dataNames=map.ESepiGen4c1l[["GeneSymbols"]])
pd.ESepiGen4c1l<-data.frame(Condition=sapply(colnames(p.ESepiGen4c1l$mRNA.Seq),
function(x) unlist(strsplit(x,'_'))[1]),stringsAsFactors=FALSE)
pd.ESepiGen4c1l$color<-c(rep("red",2),rep("green",3),rep("blue",2),rep("black",2))
names(pd.ESepiGen4c1l$color)<-pd.ESepiGen4c1l$Cond
dPCA2ESepi<- data.frame(cbind(t(PCA2ESepi[[1]]),pd.ESepiGen4c1l))
#plot pca
library(ggplot2)
setEpiCOL <- scale_colour_manual(values = c("red","green","blue","black"),
guide = guide_legend(title="Lineage"))
pPC2ESepiGen4c1l <- ggplot(dPCA2ESepi, aes(x=PC1, y=PC2, colour=Condition)) +
geom_point(size=5) + setEpiCOL +
theme(legend.position=c(0,0), legend.justification=c(0,0),
panel.background = element_rect(fill = "white"),
legend.direction = "horizontal",
plot.title = element_text(vjust = 0,hjust=0,face="bold")) +
labs(title = "Encode RNAseq in target PC1 & PC2",
x=paste("Projected PC1 (",round(PCA2ESepi[[2]][1],2),"% of varience)",sep=""),
y=paste("Projected PC2 (",round(PCA2ESepi[[2]][2],2),"% of varience)",sep=""))
library(gridExtra)
grid.arrange(pPCA,pPC2ESepiGen4c1l,nrow=1)
NMF projection
NMF decomposes a data matrix of $D$ with $N$ genes as rows and $M$ samples as columns, into two matrices, as $D ~ AP$. The pattern matrix P has rows associated with BPs in samples and the amplitude matrix A has columns indicating the relative association of a given gene, where the total number of BPs (k) is an input parameter. CoGAPS and GWCoGAPS seek a pattern matrix (${\bf{P}}$) and the corresponding distribution matrix of weights (${\bf{A}}$) whose product forms a mock data matrix (${\bf{M}}$) that represents the gene-wise data ${\bf{D}}$ within noise limits ($\boldsymbol{\varepsilon}$). That is,
\begin{equation}
{\bf{D}} = {\bf{M}} + \boldsymbol{\varepsilon} = {\bf{A}}{\bf{P}} + \boldsymbol{\varepsilon}.
\label{eq:matrixDecomp}
\end{equation}
The number of rows in ${\bf{P}}$ (columns in ${\bf{A}}$) defines the number of biological patterns (k) that CoGAPS/GWCoGAPS will infer from the number of nonorthogonal basis vectors required to span the data space. As in the Bayesian Decomposition algorithm @Ochs2006, the matrices ${\bf{A}}$ and ${\bf{P}}$ in CoGAPS are assumed to have the atomic prior described in @Sibisi1997. In the CoGAPS/GWCoGAPS implementation, $\alpha_{A}$ and $\alpha_{P}$ are corresponding parameters for the expected number of atoms which map to each matrix element in ${\bf{A}}$ and ${\bf{P}}$, respectively. The corresponding matrices ${\bf{A}}$ and ${\bf{P}}$ are found by MCMC sampling.
Projection of CoGAPS/GWCoGAPS patterns is implemented by solving the factorization in \ref{eq:matrixDecomp} for the new data matrix where ${\bf{A}}$ is the fixed nonorthogonal basis vectors comprising the average of the posterior mean for the CoGAPS/GWCoGAPS simulations performed on the original data. The patterns ${\bf{P}}$ in the new data associated with this amplitude matrix is estimated using the least-squares fit to the new data implemented with the lmFit function in the limma package (cite Limma). The \texttt{projectR} function has S4 method for class \texttt{Linear Embedding Matrix, LME}.
library(projectR)
projectR(data, loadings,dataNames = NULL, loadingsNames = NULL,
NP = NA, full = FALSE)
Input Arguments
The inputs that must be set each time are only the data and patterns, with all other inputs having default values. However, inconguities between gene names--rownames of the loadings object and either rownames of the data object will throw errors and, subsequently, should be checked before running.
The arguments are as follows:
\begin{description}
\item[data]{a target dataset to be projected into the pattern space}
\item[loadings]{a CogapsResult object}
\item[dataNames]{rownames (eg. gene names) of the target dataset, if different from existing rownames of data}
\item[loadingsNames] loadingsNames rownames (eg. gene names) of the loadings to be matched with dataNames
\item[NP]{vector of integers indicating which columns of loadings object to use. The default of NP = NA will use entire matrix.}
\item[full]{logical indicating whether to return the full model solution. By default only the new pattern object is returned.}
\end{description}
Output
The basic output of the base projectR function, i.e. \texttt{full=FALSE}, returns \texttt{projectionPatterns} representing relative weights for the samples from the new data in this previously defined feature space, or set of feature spaces. The full output of the base projectR function, i.e. \texttt{full=TRUE}, returns \texttt{projectionFit}, a list containing \texttt{projectionPatterns} and \texttt{Projection}. The \texttt{Projection} object contains additional information from the proceedure used to obtain the \texttt{projectionPatterns}. For the the the base projectR function, \texttt{Projection} is the full lmFit model from the package r BiocStyle::Biocpkg('limma').
Obtaining CoGAPS patterns to project.
# get data
library(projectR)
AP <- get(data("AP.RNAseq6l3c3t")) #CoGAPS run data
AP <- AP$Amean
# heatmap of gene weights for CoGAPs patterns
library(gplots)
pNMF<-heatmap.2(as.matrix(AP),col=bluered, trace='none',
distfun=function(c) as.dist(1-cor(t(c))) ,
cexCol=1,cexRow=.5,scale = "row",
hclustfun=function(x) hclust(x, method="average")
)
Projecting CoGAPS objects
# data to project into PCs from RNAseq6l3c3t expression data
library(projectR)
data('p.ESepiGen4c1l4')
data('p.RNAseq6l3c3t')
NMF2ESepi <- projectR(p.ESepiGen4c1l$mRNA.Seq,loadings=AP,full=TRUE,
dataNames=map.ESepiGen4c1l[["GeneSymbols"]])
dNMF2ESepi<- data.frame(cbind(t(NMF2ESepi),pd.ESepiGen4c1l))
#plot pca
library(ggplot2)
setEpiCOL <- scale_colour_manual(values = c("red","green","blue","black"),
guide = guide_legend(title="Lineage"))
pNMF2ESepiGen4c1l <- ggplot(dNMF2ESepi, aes(x=X1, y=X2, colour=Condition)) +
geom_point(size=5) + setEpiCOL +
theme(legend.position=c(0,0), legend.justification=c(0,0),
panel.background = element_rect(fill = "white"),
legend.direction = "horizontal",
plot.title = element_text(vjust = 0,hjust=0,face="bold"))
labs(title = "Encode RNAseq in target PC1 & PC2",
x=paste("Projected PC1 (",round(PCA2ESepi[[2]][1],2),"% of varience)",sep=""),
y=paste("Projected PC2 (",round(PCA2ESepi[[2]][2],2),"% of varience)",sep=""))
Clustering projection
As cannonical projection is not defined for clustering objects, the projectR package offers two transfer learning inspired methods to achieve the "projection" of clustering objects. These methods are defined by the function used to quantify and transfer the relationships which define each cluster in the original data set to the new dataset. Briefly, \texttt{cluster2pattern} uses the corelation of each genes expression to the mean of each cluster to define continuous weights. These weights are output as a pclust object which can serve as input to \texttt{projectR}. Alternatively, the \texttt{intersectoR} function can be used to test for significant overlap between two clustering objects. Both \texttt{cluster2pattern} and \texttt{intersectoR} methods are coded for a generic list structure with additional S4 class methods for kmeans and hclust objects. Further details and examples are provided in the followin respecitive sections.
cluster2pattern
\texttt{cluster2pattern} uses the corelation of each genes expression to the mean of each cluster to define continuous weights.
library(projectR)
cluster2pattern(clusters = NA, NP = NA, Data = NA)
Input Arguments
The inputs that must be set each time are the clusters and data.
The arguments are as follows:
\begin{description}
\item[clusters]{a clustering object}
\item[NP]{either the number of clusters desired or the subset of clusters to use}
\item[data]{data used to make clusters object}
\end{description}
Output
The output of the \texttt{cluster2pattern} function is a \texttt{pclust} class object; specifically, a matrix of genes (rows) by clusters (columns). A gene's value outside of its assigned cluster is zero. For the cluster containing a given gene, the gene's value is the corelation of the gene's expression to the mean of that cluster.
intersectoR
\texttt{intersectoR} function can be used to test for significant overlap between two clustering objects. The base function finds and tests the intersecting values of two sets of lists, presumably the genes associated with patterns in two different datasets. S4 class methods for hclust and kmeans objects are also available.
library(projectR)
intersectoR(pSet1 = NA, pSet2 = NA, pval = 0.05, full = FALSE, k = NULL)
Input Arguments
The inputs that must be set each time are the clusters and data.
The arguments are as follows:
\begin{description}
\item[pSet1]{a list for a set of patterns where each entry is a set of genes associated with a single pattern}
\item[pSet2]{a list for a second set of patterns where each entry is a set of genes associated with a single pattern}
\item[pval]{the maximum p-value considered significant}
\item[full]{logical indicating whether to return full data frame of signigicantly overlapping sets. Default is false will return summary matrix.}
\item[k]{numeric giving cut height for hclust objects, if vector arguments will be applied to pSet1 and pSet2 in that order}
\end{description}
Output
The output of the \texttt{intersectoR} function is a summary matrix showing the sets with stastically significant overlap under the specified p-val threshold based on a hypergeometric test. If \texttt{full==TRUE} the full data frame of signigicantly overlapping sets will also be returned.
Correlation based projection
Correlation based projection requires a matrix of gene-wise correlation values to serve as the Pattern input to the \texttt{projectR} function. This matrix can be user-generated or the result of the \texttt{correlateR} function included in the projectR package. User-generated matrixes with each row corresponding to an individual gene can be input to the generic \texttt{projectR} function. The \texttt{correlateR} function allows users to create a weight matrix for projection with values quantifying the within dataset correlation of each genes expression to the expression pattern of a particular gene or set of genes as follows.
correlateR
library(projectR)
correlateR(genes = NA, dat = NA, threshtype = "R", threshold = 0.7, absR = FALSE, ...)
Input Arguments
The inputs that must be set each time are only the genes and data, with all other inputs having default values.
The arguments are as follows:
\begin{description}
\item[genes]{gene or character vector of genes for reference expression pattern dat}
\item[data]{matrix or data frame with genes to be used for to calculate correlation}
\item[threshtype]{Default "R" indicates thresholding by R value or equivalent. Alternatively, "N" indicates a numerical cut off.}
\item[threshold]{numeric indicating value at which to make threshold}
\item[absR]{logical indicating where to include both positive and negatively correlated genes}
\item[...]{addtion imputes to the cor function}
\end{description}
Output
The output of the \texttt{correlateR} function is a \texttt{correlateR} class object. Specifically, a matrix of correlation values for those genes whose expression pattern pattern in the dataset is correlated (and anti-correlated if absR=TRUE) above the value given in as the threshold arguement. As this information may be useful in its own right, it is recommended that users inspect the \texttt{correlateR} object before using it as input to the \texttt{projectR} function.
Obtaining and visualizing \texttt{correlateR} objects.
# data to
library(projectR)
data("p.RNAseq6l3c3t")
# get genes correlated to T
cor2T<-correlateR(genes="T", dat=p.RNAseq6l3c3t, threshtype="N", threshold=10, absR=TRUE)
cor2T <- cor2T@corM
### heatmap of genes more correlated to T
indx<-unlist(sapply(cor2T,rownames))
indx <- as.vector(indx)
colnames(p.RNAseq6l3c3t)<-pd.RNAseq6l3c3t$sampleX
library(reshape2)
pm.RNAseq6l3c3t<-melt(cbind(p.RNAseq6l3c3t[indx,],indx))
library(gplots)
library(ggplot2)
library(viridis)
pCorT<-ggplot(pm.RNAseq6l3c3t, aes(variable, indx, fill = value)) +
geom_tile(colour="gray20", size=1.5, stat="identity") +
scale_fill_viridis(option="B") +
xlab("") + ylab("") +
scale_y_discrete(limits=indx) +
ggtitle("Ten genes most highly pos & neg correlated with T") +
theme(
panel.background = element_rect(fill="gray20"),
panel.border = element_rect(fill=NA,color="gray20", size=0.5, linetype="solid"),
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
axis.line = element_blank(),
axis.ticks = element_blank(),
axis.text = element_text(size=rel(1),hjust=1),
axis.text.x = element_text(angle = 90,vjust=.5),
legend.text = element_text(color="white", size=rel(1)),
legend.background = element_rect(fill="gray20"),
legend.position = "bottom",
legend.title=element_blank()
)
pCorT
Projecting correlateR objects.
# data to project into from RNAseq6l3c3t expression data
data(p.ESepiGen4c1l)
library(projectR)
cor2ESepi <- projectR(p.ESepiGen4c1l$mRNA.Seq,loadings=cor2T[[1]],full=FALSE,
dataNames=map.ESepiGen4c1l$GeneSymbols)
Try the projectR package in your browser
Any scripts or data that you put into this service are public.
projectR documentation built on Nov. 8, 2020, 7:50 p.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
knitr::opts_chunk$set(echo = TRUE) options(scipen = 1, digits = 2) set.seed(1234)
options(tinytex.verbose = TRUE)
Introduction
Technological advances continue to spur the exponential growth of biological data as illustrated by the rise of the omics—genomics, transcriptomics, epigenomics, proteomics, etc.—each with there own high throughput technologies. In order to leverage the full power of these resources, methods to integrate multiple data sets and data types must be developed. The reciprocal nature of the genomic, transcriptomic, epigenomic, and proteomic biology requires that the data provides a complementary view of cellular function and regulatory organization; however, the technical heterogeneity and massive size of high-throughput data even within a particular omic makes integrated analysis challenging. To address these challenges, we developed projectR, an R package for integrated analysis of high dimensional omic data. projectR uses the relationships defined within a given high dimensional data set, to interrogate related biological phenomena in an entirely new data set. By relying on relative comparisons within data type, projectR is able to circumvent many issues arising from technological variation. For a more extensive example of how the tools in the projectR package can be used for in silico experiments, or additional information on the algorithm, see Stein-O'Brien, et al.
Getting started with projectR
Installation Instructions
For automatic Bioconductor package installation, start R, and run:
BiocManager::install("projectR")
Methods
Projection can roughly be defined as a mapping or transformation of points from one space to another often lower dimensional space. Mathematically, this can described as a function $\varphi(x)=y : \Re^{D} \mapsto \Re^{d}$ s.t. $d \leq D$ for $x \in \Re^{D}, y \in \Re^{d}$ @Barbakh:2009bw . The projectR package uses projection functions defined in a training dataset to interrogate related biological phenomena in an entirely new data set. These functions can be the product of any one of several methods common to "omic" analyses including regression, PCA, NMF, clustering. Individual sections focussing on one specific method are included in the vignette. However, the general design of the projectR function is the same regardless.
The base projectR function
The generic projectR function is executed as follows:
library(projectR) projectR(data, loadings, dataNames=NULL, loadingsNames=NULL, NP = NULL, full = false)
Input Arguments
The inputs that must be set each time are only the data and loadings, with all other inputs having default values. However, incongruities in the feature mapping between the data and loadings, i.e. a different format for the rownames of each object, will throw errors or result in an empty mapping and should be checked before running. To overcoming mismatched feature names in the objects themselve, the /code{dataNames} and /code{loadingNames} arguements can be manually supplied by the user.
The arguments are as follows: \begin{description} \item[data]{a dataset to be projected into the pattern space} \item[loadings]{a matrix of continous values with unique rownames to be projected} \item[dataNames]{a vector containing unique name, i.e. gene names, for the rows of the target dataset to be used to match features with the loadings, if not provided by \texttt{rownames(data)}. Order of names in vector must match order of rows in data.} \item[loadingsNames]{a vector containing unique names, i.e. gene names, for the rows of loadings to be used to match features with the data, if not provided by \texttt{rownames(loadings)}. Order of names in vector must match order of rows in loadings.} \item[NP]{vector of integers indicating which columns of loadings object to use. The default of NP = NA will use entire matrix.} \item[full]{logical indicating whether to return the full model solution. By default only the new pattern object is returned.} \end{description}
The \texttt{loadings} arguement in the generic projectR function is suitable for use with any genernal feature space, or set of feature spaces, whose rows annotation links them to the data to be projected. Ex: the coeffients associated with individual genes as the result of regression analysis or the amplituded values of individual genes as the result of non-negative matrix factorization (NMF).
Output
The basic output of the base projectR function, i.e. \texttt{full=FALSE}, returns \texttt{projectionPatterns} representing relative weights for the samples from the new data in this previously defined feature space, or set of feature spaces. The full output of the base projectR function, i.e. \texttt{full=TRUE}, returns \texttt{projectionFit}, a list containing \texttt{projectionPatterns} and \texttt{Projection}. The \texttt{Projection} object contains additional information from the proceedure used to obtain the \texttt{projectionPatterns}. For the the the base projectR function, \texttt{Projection} is the full lmFit model from the package r BiocStyle::Biocpkg("limma")
PCA projection
Projection of principle components is achieved by matrix multiplication of a new data set by previously generated eigenvectors, or gene loadings. If the original data were standardized such that each gene is centered to zero average expression level, the principal components are normalized eigenvectors of the covariance matrix of the genes. Each PC is ordered according to how much of the variation present in the data they contain. Projection of the original samples into each PC will maximized the variance of the samples in the direction of that component and uncorrelated to previous components. Projection of new data places the new samples into the PCs defined by the original data. Because the components define an orthonormal basis set, they provide an isomorphism between a vector space, $V$, and $\Re^n$ which preserves inner products. If $V$ is an inner product space over $\Re$ with orthonormal basis $B = v_1,...,v_n$ and $v \epsilon V s.t [v]_B = (r_1,...,r_n)$, then finding the coordinate of $v_i$ in $v$ is precisely the inner product of $v$ with $v_i$, i.e. $r_i = \langle v,v_i \rangle$. This formulation is implemented for only those genes belonging to both the new data and the PC space. The \texttt{projectR} function has S4 method for class \texttt{prcomp}.
Obtaining PCs to project.
# data to define PCs library(projectR) data(p.RNAseq6l3c3t) # do PCA on RNAseq6l3c3t expression data pc.RNAseq6l3c3t<-prcomp(t(p.RNAseq6l3c3t)) pcVAR <- round(((pc.RNAseq6l3c3t$sdev)^2/sum(pc.RNAseq6l3c3t$sdev^2))*100,2) dPCA <- data.frame(cbind(pc.RNAseq6l3c3t$x,pd.RNAseq6l3c3t)) #plot pca library(ggplot2) setCOL <- scale_colour_manual(values = c("blue","black","red"), name="Condition:") setFILL <- scale_fill_manual(values = c("blue","black","red"),guide = FALSE) setPCH <- scale_shape_manual(values=c(23,22,25,25,21,24),name="Cell Line:") pPCA <- ggplot(dPCA, aes(x=PC1, y=PC2, colour=ID.cond, shape=ID.line, fill=ID.cond)) + geom_point(aes(size=days),alpha=.6)+ setCOL + setPCH + setFILL + scale_size_area(breaks = c(2,4,6), name="Day") + theme(legend.position=c(0,0), legend.justification=c(0,0), legend.direction = "horizontal", panel.background = element_rect(fill = "white",colour=NA), legend.background = element_rect(fill = "transparent",colour=NA), plot.title = element_text(vjust = 0,hjust=0,face="bold")) + labs(title = "PCA of hPSC PolyA RNAseq", x=paste("PC1 (",pcVAR[1],"% of varience)",sep=""), y=paste("PC2 (",pcVAR[2],"% of varience)",sep=""))
Projecting prcomp objects
# data to project into PCs from RNAseq6 l3c3t expression data data(p.ESepiGen4c1l) library(projectR) PCA2ESepi <- projectR(data = p.ESepiGen4c1l$mRNA.Seq,loadings=pc.RNAseq6l3c3t, full=TRUE, dataNames=map.ESepiGen4c1l[["GeneSymbols"]]) pd.ESepiGen4c1l<-data.frame(Condition=sapply(colnames(p.ESepiGen4c1l$mRNA.Seq), function(x) unlist(strsplit(x,'_'))[1]),stringsAsFactors=FALSE) pd.ESepiGen4c1l$color<-c(rep("red",2),rep("green",3),rep("blue",2),rep("black",2)) names(pd.ESepiGen4c1l$color)<-pd.ESepiGen4c1l$Cond dPCA2ESepi<- data.frame(cbind(t(PCA2ESepi[[1]]),pd.ESepiGen4c1l)) #plot pca library(ggplot2) setEpiCOL <- scale_colour_manual(values = c("red","green","blue","black"), guide = guide_legend(title="Lineage")) pPC2ESepiGen4c1l <- ggplot(dPCA2ESepi, aes(x=PC1, y=PC2, colour=Condition)) + geom_point(size=5) + setEpiCOL + theme(legend.position=c(0,0), legend.justification=c(0,0), panel.background = element_rect(fill = "white"), legend.direction = "horizontal", plot.title = element_text(vjust = 0,hjust=0,face="bold")) + labs(title = "Encode RNAseq in target PC1 & PC2", x=paste("Projected PC1 (",round(PCA2ESepi[[2]][1],2),"% of varience)",sep=""), y=paste("Projected PC2 (",round(PCA2ESepi[[2]][2],2),"% of varience)",sep=""))
library(gridExtra) grid.arrange(pPCA,pPC2ESepiGen4c1l,nrow=1)
NMF projection
NMF decomposes a data matrix of $D$ with $N$ genes as rows and $M$ samples as columns, into two matrices, as $D ~ AP$. The pattern matrix P has rows associated with BPs in samples and the amplitude matrix A has columns indicating the relative association of a given gene, where the total number of BPs (k) is an input parameter. CoGAPS and GWCoGAPS seek a pattern matrix (${\bf{P}}$) and the corresponding distribution matrix of weights (${\bf{A}}$) whose product forms a mock data matrix (${\bf{M}}$) that represents the gene-wise data ${\bf{D}}$ within noise limits ($\boldsymbol{\varepsilon}$). That is, \begin{equation} {\bf{D}} = {\bf{M}} + \boldsymbol{\varepsilon} = {\bf{A}}{\bf{P}} + \boldsymbol{\varepsilon}. \label{eq:matrixDecomp} \end{equation} The number of rows in ${\bf{P}}$ (columns in ${\bf{A}}$) defines the number of biological patterns (k) that CoGAPS/GWCoGAPS will infer from the number of nonorthogonal basis vectors required to span the data space. As in the Bayesian Decomposition algorithm @Ochs2006, the matrices ${\bf{A}}$ and ${\bf{P}}$ in CoGAPS are assumed to have the atomic prior described in @Sibisi1997. In the CoGAPS/GWCoGAPS implementation, $\alpha_{A}$ and $\alpha_{P}$ are corresponding parameters for the expected number of atoms which map to each matrix element in ${\bf{A}}$ and ${\bf{P}}$, respectively. The corresponding matrices ${\bf{A}}$ and ${\bf{P}}$ are found by MCMC sampling.
Projection of CoGAPS/GWCoGAPS patterns is implemented by solving the factorization in \ref{eq:matrixDecomp} for the new data matrix where ${\bf{A}}$ is the fixed nonorthogonal basis vectors comprising the average of the posterior mean for the CoGAPS/GWCoGAPS simulations performed on the original data. The patterns ${\bf{P}}$ in the new data associated with this amplitude matrix is estimated using the least-squares fit to the new data implemented with the lmFit function in the limma package (cite Limma). The \texttt{projectR} function has S4 method for class \texttt{Linear Embedding Matrix, LME}.
library(projectR)
projectR(data, loadings,dataNames = NULL, loadingsNames = NULL,
NP = NA, full = FALSE)
Input Arguments
The inputs that must be set each time are only the data and patterns, with all other inputs having default values. However, inconguities between gene names--rownames of the loadings object and either rownames of the data object will throw errors and, subsequently, should be checked before running.
The arguments are as follows: \begin{description} \item[data]{a target dataset to be projected into the pattern space} \item[loadings]{a CogapsResult object} \item[dataNames]{rownames (eg. gene names) of the target dataset, if different from existing rownames of data} \item[loadingsNames] loadingsNames rownames (eg. gene names) of the loadings to be matched with dataNames \item[NP]{vector of integers indicating which columns of loadings object to use. The default of NP = NA will use entire matrix.} \item[full]{logical indicating whether to return the full model solution. By default only the new pattern object is returned.} \end{description}
Output
The basic output of the base projectR function, i.e. \texttt{full=FALSE}, returns \texttt{projectionPatterns} representing relative weights for the samples from the new data in this previously defined feature space, or set of feature spaces. The full output of the base projectR function, i.e. \texttt{full=TRUE}, returns \texttt{projectionFit}, a list containing \texttt{projectionPatterns} and \texttt{Projection}. The \texttt{Projection} object contains additional information from the proceedure used to obtain the \texttt{projectionPatterns}. For the the the base projectR function, \texttt{Projection} is the full lmFit model from the package r BiocStyle::Biocpkg('limma').
Obtaining CoGAPS patterns to project.
# get data library(projectR) AP <- get(data("AP.RNAseq6l3c3t")) #CoGAPS run data AP <- AP$Amean # heatmap of gene weights for CoGAPs patterns library(gplots) pNMF<-heatmap.2(as.matrix(AP),col=bluered, trace='none', distfun=function(c) as.dist(1-cor(t(c))) , cexCol=1,cexRow=.5,scale = "row", hclustfun=function(x) hclust(x, method="average") )
Projecting CoGAPS objects
# data to project into PCs from RNAseq6l3c3t expression data library(projectR) data('p.ESepiGen4c1l4') data('p.RNAseq6l3c3t') NMF2ESepi <- projectR(p.ESepiGen4c1l$mRNA.Seq,loadings=AP,full=TRUE, dataNames=map.ESepiGen4c1l[["GeneSymbols"]]) dNMF2ESepi<- data.frame(cbind(t(NMF2ESepi),pd.ESepiGen4c1l)) #plot pca library(ggplot2) setEpiCOL <- scale_colour_manual(values = c("red","green","blue","black"), guide = guide_legend(title="Lineage")) pNMF2ESepiGen4c1l <- ggplot(dNMF2ESepi, aes(x=X1, y=X2, colour=Condition)) + geom_point(size=5) + setEpiCOL + theme(legend.position=c(0,0), legend.justification=c(0,0), panel.background = element_rect(fill = "white"), legend.direction = "horizontal", plot.title = element_text(vjust = 0,hjust=0,face="bold")) labs(title = "Encode RNAseq in target PC1 & PC2", x=paste("Projected PC1 (",round(PCA2ESepi[[2]][1],2),"% of varience)",sep=""), y=paste("Projected PC2 (",round(PCA2ESepi[[2]][2],2),"% of varience)",sep=""))
Clustering projection
As cannonical projection is not defined for clustering objects, the projectR package offers two transfer learning inspired methods to achieve the "projection" of clustering objects. These methods are defined by the function used to quantify and transfer the relationships which define each cluster in the original data set to the new dataset. Briefly, \texttt{cluster2pattern} uses the corelation of each genes expression to the mean of each cluster to define continuous weights. These weights are output as a pclust object which can serve as input to \texttt{projectR}. Alternatively, the \texttt{intersectoR} function can be used to test for significant overlap between two clustering objects. Both \texttt{cluster2pattern} and \texttt{intersectoR} methods are coded for a generic list structure with additional S4 class methods for kmeans and hclust objects. Further details and examples are provided in the followin respecitive sections.
cluster2pattern
\texttt{cluster2pattern} uses the corelation of each genes expression to the mean of each cluster to define continuous weights.
library(projectR) cluster2pattern(clusters = NA, NP = NA, Data = NA)
Input Arguments
The inputs that must be set each time are the clusters and data.
The arguments are as follows: \begin{description} \item[clusters]{a clustering object} \item[NP]{either the number of clusters desired or the subset of clusters to use} \item[data]{data used to make clusters object} \end{description}
Output
The output of the \texttt{cluster2pattern} function is a \texttt{pclust} class object; specifically, a matrix of genes (rows) by clusters (columns). A gene's value outside of its assigned cluster is zero. For the cluster containing a given gene, the gene's value is the corelation of the gene's expression to the mean of that cluster.
intersectoR
\texttt{intersectoR} function can be used to test for significant overlap between two clustering objects. The base function finds and tests the intersecting values of two sets of lists, presumably the genes associated with patterns in two different datasets. S4 class methods for hclust and kmeans objects are also available.
library(projectR) intersectoR(pSet1 = NA, pSet2 = NA, pval = 0.05, full = FALSE, k = NULL)
Input Arguments
The inputs that must be set each time are the clusters and data.
The arguments are as follows: \begin{description} \item[pSet1]{a list for a set of patterns where each entry is a set of genes associated with a single pattern} \item[pSet2]{a list for a second set of patterns where each entry is a set of genes associated with a single pattern} \item[pval]{the maximum p-value considered significant} \item[full]{logical indicating whether to return full data frame of signigicantly overlapping sets. Default is false will return summary matrix.} \item[k]{numeric giving cut height for hclust objects, if vector arguments will be applied to pSet1 and pSet2 in that order} \end{description}
Output
The output of the \texttt{intersectoR} function is a summary matrix showing the sets with stastically significant overlap under the specified p-val threshold based on a hypergeometric test. If \texttt{full==TRUE} the full data frame of signigicantly overlapping sets will also be returned.
Correlation based projection
Correlation based projection requires a matrix of gene-wise correlation values to serve as the Pattern input to the \texttt{projectR} function. This matrix can be user-generated or the result of the \texttt{correlateR} function included in the projectR package. User-generated matrixes with each row corresponding to an individual gene can be input to the generic \texttt{projectR} function. The \texttt{correlateR} function allows users to create a weight matrix for projection with values quantifying the within dataset correlation of each genes expression to the expression pattern of a particular gene or set of genes as follows.
correlateR
library(projectR) correlateR(genes = NA, dat = NA, threshtype = "R", threshold = 0.7, absR = FALSE, ...)
Input Arguments
The inputs that must be set each time are only the genes and data, with all other inputs having default values.
The arguments are as follows: \begin{description} \item[genes]{gene or character vector of genes for reference expression pattern dat} \item[data]{matrix or data frame with genes to be used for to calculate correlation} \item[threshtype]{Default "R" indicates thresholding by R value or equivalent. Alternatively, "N" indicates a numerical cut off.} \item[threshold]{numeric indicating value at which to make threshold} \item[absR]{logical indicating where to include both positive and negatively correlated genes} \item[...]{addtion imputes to the cor function} \end{description}
Output
The output of the \texttt{correlateR} function is a \texttt{correlateR} class object. Specifically, a matrix of correlation values for those genes whose expression pattern pattern in the dataset is correlated (and anti-correlated if absR=TRUE) above the value given in as the threshold arguement. As this information may be useful in its own right, it is recommended that users inspect the \texttt{correlateR} object before using it as input to the \texttt{projectR} function.
Obtaining and visualizing \texttt{correlateR} objects.
# data to library(projectR) data("p.RNAseq6l3c3t") # get genes correlated to T cor2T<-correlateR(genes="T", dat=p.RNAseq6l3c3t, threshtype="N", threshold=10, absR=TRUE) cor2T <- cor2T@corM ### heatmap of genes more correlated to T indx<-unlist(sapply(cor2T,rownames)) indx <- as.vector(indx) colnames(p.RNAseq6l3c3t)<-pd.RNAseq6l3c3t$sampleX library(reshape2) pm.RNAseq6l3c3t<-melt(cbind(p.RNAseq6l3c3t[indx,],indx)) library(gplots) library(ggplot2) library(viridis) pCorT<-ggplot(pm.RNAseq6l3c3t, aes(variable, indx, fill = value)) + geom_tile(colour="gray20", size=1.5, stat="identity") + scale_fill_viridis(option="B") + xlab("") + ylab("") + scale_y_discrete(limits=indx) + ggtitle("Ten genes most highly pos & neg correlated with T") + theme( panel.background = element_rect(fill="gray20"), panel.border = element_rect(fill=NA,color="gray20", size=0.5, linetype="solid"), panel.grid.major = element_blank(), panel.grid.minor = element_blank(), axis.line = element_blank(), axis.ticks = element_blank(), axis.text = element_text(size=rel(1),hjust=1), axis.text.x = element_text(angle = 90,vjust=.5), legend.text = element_text(color="white", size=rel(1)), legend.background = element_rect(fill="gray20"), legend.position = "bottom", legend.title=element_blank() )
pCorT
Projecting correlateR objects.
# data to project into from RNAseq6l3c3t expression data data(p.ESepiGen4c1l) library(projectR) cor2ESepi <- projectR(p.ESepiGen4c1l$mRNA.Seq,loadings=cor2T[[1]],full=FALSE, dataNames=map.ESepiGen4c1l$GeneSymbols)
Try the projectR package in your browser
Any scripts or data that you put into this service are public.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.