BACKGROUND
there has been rapid progress in high-throughput technologies and platforms to assay global mrna, mirna, methylation, proteins, and metabolite profiles of cells are readily available.
advances in rna-sequencing and mass spectrometry  based proteomics have dramatically improved coverage and quality of genomic, transcriptomic and proteomic profiling  <cit> .
increasing number of studies including the cancer genome atlas  and encyclopedia of dna elements  projects systematically profile large number of biological samples resulting in multiple levels of quantitative information  <cit> .
recent advances of ms based proteomics provide a complementary approach to genomics and transcriptomic technologies  <cit>  and systematic analyses can now be carried out to identify and quantify the majority of proteins expressed in human cells  <cit> .
these data yield unprecedented views of molecular building blocks and the machinery of cells.
interpreting these large-scale datasets to derive information about a biological system represents a considerable challenge often faced by investigators.
multiple omics data analysis can be broadly defined by some common questions, which are dependent on the data collected; multiple datasets measuring the same biological molecules or multiple datasets each measuring different biological molecules.
in the first case, given multiple transcriptomics data from different microarray or rna-sequencing studies, the aim may be to discover which platform is the most informative with highest quality data, identify robust biomarkers across datasets or highlight platform specific discrepancies in measurements.
in the second case, given multiple different data such as transcripts, proteins and metabolites, the objective may be to integrate and concatenate information to increase the breath and coverage of available data in a biological network.
in this case, specific platform discrepancies are less important and performance of data integration is more likely to be assessed using system biology or pathway approaches.
nevertheless, both analyses face common challenges associated with integrating data from disparate technologies.
several meta-analysis studies map identifiers from each platform to a common set of identifiers to generate a single concatenated matrix for subsequent analysis  <cit> .
however, this data simplification overlooks several fundamental platform and biological biases.
platforms are not universal and measure different molecules.
filtering genes to their intersection may considerably reduce data coverage.
in addition, the many-to-many mapping of gene identifiers from multiple platforms complicates direct comparison of molecules across multiple levels.
moreover, because correlations between different platforms are probably lower than expected  <cit> , it may not provide gains in data quality or study power.
such filtering may also introduce bias because platform discrepancies could reflect biological variation.
for instance, poor correlation between a transcript and its translated protein may result from biological processes such as microrna post-transcriptional repression  <cit> .
similarly, correlations between proteins and metabolites of pathways can diverge if proteins are expressed in an inactive form, in which case its abundance may not represent activity.
ordination methods, such as principal component analysis , independent component analysis  and correspondence analysis , are exploratory data analysis approaches that have been applied to analyze omics data including transcriptome and proteome studies  <cit> .
graphical representation of measurements  and variables  on a lower dimensional space facilitates interpretation of global variance structure and identification of the most informative  features across datasets.
these methods permit visualization of data that have considerable levels of noise and data where the number of variables exceeds the number of measurements, which is typical in omics studies.
however, these approaches do not solve the problem of comparing many datasets simultaneously.
studies have extended these approaches to couple two datasets together  <cit> .
one such approach is co-inertia analysis   <cit> .
cia was originally applied to study ecological and environmental tables, where it was employed to link environmental variables with species characteristics  <cit> .
culhane and colleagues introduced cia in genomics, when they compared data from two microarray platforms  <cit> .
an advantage of this method is that it does not require the mapping or filtering of genes to a common set.
the relationship between co-inertia analysis and related methods including procrustes analysis  <cit> , canonical correlation analysis with elastic net penalization  and sparse partial least squares  have been described previously  <cit> .
cia and spls both maximize the covariance between eigenvectors and are efficient in determining main individual effects in paired dataset analysis.
by contrast cca-en maximizes the correlation between eigenvectors and tends to discover effects present in both datasets, but may omit to discover strong individual effects.
variables selected by cca-en and spls are highly similar but cia selected marginally different marker genes that may have some redundancy  <cit> .
a noteworthy advantage of cia is that it can be coupled with several dimension reduction approaches, including pca or correspondence analysis, such that it can accommodate both discrete count data  and continuous data  <cit> .
these approaches are performed on each dataset separately and can be integrated using cia  <cit> .
however, all above methods including cia are limited to the analysis of two datasets, limiting their application in modern multi-omics studies.
several approaches have been proposed for integrating more than two datasets, such as consensus pca   <cit> , regularized generalized canonical correlation analysis   <cit> , sparse generalized canonical correlation analysis   <cit>  and penalized canonical correlation analysis   <cit> .
sgcca and pcca originally focus on the feature selection from multiple datasets, but also can be used for multiple table integration problem.
here, we describe another method, multiple co-inertia analysis , for the analysis of more than two omics datasets, extending its application in the field of environmental science and, recently, phylogenetics  <cit> .
mcia is related to consensus pca  which both maximize the square covariance between eigenvectors and are subject to similar constraints  <cit> .
cpca is less sensitive to multi-collinearity within each dataset than generalized canonical correlation analysis  <cit> .
we illustrate the application of mcia using two different examples, and show that integrated analysis is more insightful than analysis of the individual datasets.
first, we demonstrated the power of mcia via applying it to the integration and comparison of multi-omics data independent of data annotation.
we employed mcia to identify common relationships among multiple gene and protein expression data of the nci- <dig> cancer cell line panel of the national cancer institute  <cit> .
the integrated analysis revealed that cell lines are clustered according to anatomical tissue source and showed a significant degree of correlation between transcript and protein expression.
second, we assessed the concordance in gene expression data obtained from microarray and next generation rna-sequencing of  <dig> samples of high grade serous ovarian cancer.
mcia integrated ovarian cancer gene expression data from different sources which captured distinct subsets of the transcriptome  to reveal a set of biomarkers that were consistently highly ranked by all four platforms and were biologically relevant to ovarian cancer.
to enable community access to mcia, we implemented the method into the r-bioconductor  package as an easy-to-use tool for bioinformaticians and biologists.
methods
mathematical basis of mcia
a typical omics dataset is a matrix where the number of features exceeds the number of measurements .
mcia requires a set of tables where either features or measurements are matched and have equal weights.
mcia is performed in a two-step process.
first a one table ordination method, such as pca, coa or non-symmetric correspondence analysis   <cit>  is applied on each dataset separately, which transforms data into comparable lower dimensional spaces.
in our analysis, given an omics data table m =  with 1 ≤ i ≤ n and 1 ≤ j ≤ q, where m is a  matrix, i indicates row index and j for column index.
we denote the row and column sums of m as mi+ and m+j respectively, and m++ as the grand total.
the relative contribution or weight of row i to the total variation in the data set is denoted ri and calculated as ri = mi+/m++ while the relative contribution of column j is denoted as cj = m+j/m++.
similarly, the contribution of each individual element of m to the total variation pij can be calculated as pij = mij/m++.
we then derive a new matrix x with the values defined above as
xij=pijri−cj
where xij is the centered row profile, i.e.
the relative abundance of selected variable to the measurement’s weight.
the second step in mcia is a generalization of cia  <cit> .
it solves the problem of simultaneous analysis of a set of statistical triplet  where k =  <dig>  …, k,…, k and xk is a set of transformed matrices.
qk is a qk × qk matrix with rij in diagonal elements, indicating the hyperspace of features metrics.
d is an n × n matrix which is an identity matrix indicating equal weight across all columns in all tables.
mcia maximizes the sum of the squared covariance between scores of each table with synthetic axes ν, that is:
fu <dig> …,uk,…,uk,v=∑k=1kwkcov2xkqkuk,v
where cov <dig> stands for the square of covariance of quantities inside parenthesis and ωk is the weight of each table.
the v represents the reference structure or synthetic center and uk are auxiliary axes.
the score of each individual table would then be xkqkuk.
in contrast with other ordination methods, mcia finds solutions  sequentially.
multiple matrices xk can be weighted and concatenated to a single matrix x = .
similarly, a single feature metric q could be concatenated as q = .the first order solutions of u <dig> to uk <dig> and v <dig> are given by the first principal component of the following eigen-system:
wxqxtdv=λv
then the normalized auxiliary axis uk <dig> are
uk1=xktdv1/||xktdv1||qkk= <dig> …,k
where ||•|| is the norm in the qk metric.
the subsequent solutions are found with residual matrices from the calculation of the first order solution with the constraint that the remaining order axes are orthogonal with the previous sets, namely:
vjtdvs=0andukjtqkuks=01≤j<s
the residual matrices used by second order solution is deflated as
x1order2=x1‒x1pk <dig>
where the projection matrix pk <dig> is
pk1=uk1uk1qkuk1t‒1uk1qk
the superscript t and − <dig> stand for matrix transposition and matrix inversion respectively.
therefore, the formula  removes the dimension that is spanned by vector uk <dig>  to get a residual matrix, which is passed to the svd to find the second order solution.
these steps are repeated until the desired number of axes  is generated.
as a result, mcia provides a simultaneous ordination of columns  and rows  of multiple tables within the same hyperspace, with features or measurements sharing similar trends will be closely projected.
the detailed description of mcia and the proof that these axes are maximally co-variant are given in chessel and hanafi  <cit> .
datasets
we analyzed publicly available sets of data from two studies:  transcriptomic  <cit>  and proteomic  <cit>  datasets of the nci- <dig> cancer cell line panel, the latter one generated in our group, and  an ovarian cancer dataset generated as part of the tcga project  <cit> .
in each study, there are multiple datasets measuring molecules  from the same samples .
nci- <dig> data
the nci- <dig> panel is a collection of  <dig> cancer cell lines of leukemia, lymphomas, melanomas and carcinomas of ovarian, renal, breast, prostate, colon, lung and central nervous system  origin.
the nci- <dig> transcriptome data were downloaded from cellminer  <cit>  and were obtained on four different platforms; affymetrix hg-u <dig> plus  <dig> , hg-u <dig>  hg-u <dig> and agilent ge 4x44k  <cit> .
affymetrix data were normalized using gc robust multichip averaging gcrma;  <cit>  and agilent data were log transformed as obtained from the cellminer.
although data filtering is not required to perform mcia, to facilitate data interpretation, microarray data were filtered to exclude probes that do not map to an official hugo gene symbol.
the probe with highest average value was retained when multiple probes mapped to the same gene.
filtering produced datasets of  <dig> ;  <dig> ;  <dig>  and  <dig>  genes on agilent, hg-u <dig>  hg-u <dig> and hg-u <dig> plus  <dig>  platforms respectively.
the lung cancer cell line nci-h <dig> was excluded since its expression profile was not available on the hg-u <dig> platform.
a venn diagram representing the overlapping genes in the processed data for each platform is provided in additional file 1: figure s <dig>
the proteome profiles of cell lines were produced from a conventional gelc-ms/ms approach and label-free quantification, as described in  <cit> .
the international protein index  identifiers were mapped to official gene symbol to facilitate subsequent pathway interpretation.
data were log transformed  and no filtering or additional normalization were performed.
this dataset represents  <dig>  protein expressions across  <dig> cell line in nci- <dig> panel.
ovarian cancer datasets
gene expression of tumors from ovarian cancer patients were profiled using two microarray platforms  and rna-sequencing on illumina hiseq platform.
data were downloaded from the nci-tcga data portal 07/08/2013;  <cit> .
patient samples  that were present in all four datasets were included in the analysis.
the agilent and affymetrix data were normalized and summarized by lowess and robust multichip averaging , respectively  <cit> .
the transcript expression levels of the illumina rna-sequencing data were determined using two different pre-processing pipelines  denoted as rnaseq and rnaseqv <dig>  respectively.
normalization and quantification of rnaseq followed the rpkm method  <cit>  whereas the alignment and gene expression quantification in rnaseqv <dig> were obtained by mapsplice and rsem  <cit> .
in rnaseq and rnaseqv2;  <dig>  and  <dig>  genes were detected .
these data were filtered to exclude genes with more than  <dig> missing values.
only genes mapped to an official gene symbol were retained.
for the features mapped to the same gene symbol, the one with the largest average expression value was kept.
remaining missing values were replaced with a positive value far smaller than the lowest expression value in each dataset  and then, the expression values were log transformed .
after filtering, the agilent, affymetrix, rnaseq and rnaseqv <dig> datasets contained  <dig> ;  <dig> ;  <dig>  and  <dig>  gene expression measurements respectively.
the venn diagram representing the overlap of genes in these datasets is shown in additional file 1: figure s <dig>
RESULTS
integrated analysis of the nci- <dig> cell line transcriptome and proteome
the nci- <dig> panel, a collection of  <dig> cancer cell lines derived from nine different tissues  has been extensively used in in vitro high-throughput drug screen assays.
they have been molecularly profiled using comparative genomic hybridization array  <cit> , karyotype analysis  <cit> , dna mutational analysis  <cit> , transcripts expression array  <cit> , microarrays for microrna expression  <cit>  and protein expression  <cit> .
mcia was applied as an exploratory analysis of four transcriptomic studies  and one proteomic study  of the  <dig> cell lines.
figure 1a shows the projection of cell lines onto the first two principal components  of mcia.
similar to the visualization employed in cia  <cit> , the datasets are transformed into the same projection.
the coordinates of the five measurements for each cell line are connected by lines.
the length of which indicates the divergence  between the mrnas and protein expression levels for a particular cell line.
the mcia plot of the first two principal components shows similar trends in transcriptome and proteome profiles, indicating that the most variant sources of biological information were similar.
cell lines originating from the same or closely related anatomical source of tissue were projected close to each other and converged into groups.
the colon, leukemia, melanoma, cns, renal and ovarian cell lines segregated largely according to their tissue of origin.
seven out of eight melanoma lines clustered together, and the remaining one, lox-imvi, has been reported to lack melanin production  <cit> .
these results are consistent with independently performed hierarchical clustering analysis .
there was greater divergence in the cell lines from tumors with more intrinsic molecular heterogeneity .
the transcriptome and proteome profiles of the individual breast and lung cell lines were projected close in space demonstrating that the expression profiles shared a high degree of consensus.
the tight projection of multiple data types from diverse technology platforms provides evidence that the observed spread of cell lines reflected the biological variance , as opposed to inter-study technical or other stochastic variance.
for instance, we observed that the estrogen receptor positive breast cancer cell line mcf <dig> displays an epithelial phenotype and clustered to colon cancer lines.
in contrast, the cell line negative for the estrogen receptor, hs578t, clustered with the stromal/mesenchymal cluster of glioblastoma and renal tumor cell lines.
this suggests that hs578t exhibits more invasive mesenchymal features compared to t47d and mcf <dig>
overall co-structure comparison using mcia
each pc has an associated eigenvalue which represents the amount of variability contained in that pc.
the first three pcs of the mcia accounted for  <dig> %,  <dig> % and  <dig> % of variance respectively .
the observation that the first two pcs capture less than a third of the structure in the datasets  reflects the complexity inherent in cell lines of  <dig> tumors from nine different organs.
in order to identify the contribution of each dataset to the total variance, that is, to what extent each dataset deviates or agrees with what the majority of datasets support, we extracted the mcia pseudo-eigenvalues.
figure 1b shows the pseudo-eigenvalues associated with the first two principal components of each dataset.
comparison within microarray data revealed that affymetrix hgu <dig> plus  <dig>  accounts for the highest variance on axis  <dig> and  <dig>  possibly because this platform contains informative features, or features that are poorly detected or absent on other platforms.
we observed that the similarity within transcriptome datasets is greater than the similarity between transcriptome and proteome data, which is consistent with the results shown by figure 1a.
one of the most attractive features of mcia is that it can be used to highlight lack or presence of co-structure between datasets, thus it allows selection of the strongest features from each dataset for subsequent analysis.
for instance, we observed in particular, large variation between the protein and transcript expression patterns of two cell lines, melanoma skmel <dig> and ovarian igrov <dig>  the proteome coordinates of skmel <dig> were close to the origin and far from the transcriptomic data that was projected on the negative end of pc <dig> with the other melanoma cell line data.
the poor information content in proteome data of the skmel <dig> cell line could reflect the lack of expression of melanin related genes on protein level.
similarly, the incongruence of the proteome and transcriptome data of the ovarian cell line igrov <dig> may be due to expression of less epithelial markers that projected on the positive direction of axis  <dig>
to characterize the overall correlation between each pair of high dimensional data we calculated the pair-wise rv coefficient, a multivariate generalization of the squared pearson correlation coefficient  <cit> .
for each pair of datasets, the rv-coefficient is calculated as the total co-inertia  divided by the square root of the product of the squared total inertia  from the individual analysis.
as the co-structure between two datasets increases, the rv score move towards to  <dig>  a zero rv score indicates no similarity.
the overall similarity in structure between microarray data was higher than the similarity between microarray and proteomics data; average rv coefficient >  <dig>  and  <dig>  respectively .
when mcia was performed on the same transcriptome data and the subset of proteome data that were quantified in all  <dig> cell lines , the filtered proteome data had a higher consensus to the co-structure and increased pseudo-eigenvalues .
mcia axes describe biological properties
in contrast to traditional clustering methods, mcia projects the original data onto a lower dimensional space, maximizing the covariance of each dataset with respect to the reference structure.
in mcia plots, a gene that is highly expressed in a certain cell line will be projected in the direction of this cell line and the greater the distance from the origin, the stronger the association.
in order to identify biomarkers that are highly associated with cancer cell lines of different origins, we examined the feature space of mrnas and proteins that were projected in the same direction and space .
the first axis , which explains the largest variance, separated cells with epithelial or mesenchymal characteristics, suggesting that epithelial-mesenchymal transition  is an essential mechanism defining different classes of cancers .
emt has been shown to play an important role in epithelial cell malignancy and metastasis  <cit> .
epithelial markers, including slc27a <dig>  cdh <dig>  spint <dig>  s100p and epcam had high weights on the positive side of pc <dig> .
at the opposite end, mesenchymal and collagen markers, including gja <dig>  which is involved in epicardial to mesenchymal cell transition, and tgfβ <dig> were observed .
the second  axis, pc <dig>  clearly separated melanoma and leukemia from other epithelial cancer types.
the strongest determinant of the vertical axis is a set of melanoma-related genes, namely melanoma-associated antigens , melanogenic enzyme  as well as tyr, dct, tyrp <dig>  malana, s100b and bcl2a <dig>  the top  <dig> genes with greatest weights on pc <dig> and pc <dig> were selected from each dataset  and the full list of markers is provided in additional file 2: table s <dig>  among  <dig>  selected genes,  <dig> were measured in three or more datasets.
mcia enables the study of the union of features from all studies.
among the nci- <dig> datasets, less than 12% of the total  <dig>  genes studied were measured in all five datasets.
by observing highly ranked genes across studies, one can identify robust markers that could be subject to further analysis.
integration of proteomics and transcriptomics complements the biological information
to further evaluate the biological significance of the features selected by mcia, we employed ingenuity pathway analysis  to discover significant canonical pathways which discriminate different cell lines .
in mcia plots, samples and features are projected onto the same space.
the genes with strongest association to a cell line are those projected in the same direction and have the highest weights .
as features have been transformed on the same scale, the union of features from each individual dataset can be easily extracted and concatenated to provide greater coverage in pathway analysis.
features strongly associated to each tissue type from both transcriptome and proteome datasets can be concatenated and mapped to signaling pathways.
there is no requirement to extract equal numbers of features from each data type.
for example we observed that features strongly associated with leukemia related features tended to be enriched in the proteins .
the most extreme features associated with the leukemia cell lines were selected from all platforms using their coordinates and were subjected to the functional and pathway analysis.
the full list of features, the coordinate feature selection criteria and their functional and pathway analysis are provided in additional file 3: tables s <dig> and additional file 4: table s <dig>
complementary information can be obtained by integrating data from different platforms and data types which increases the genome coverage and power of subsequent pathway analysis.
while numerous genes were over-expressed in both the transcriptome and proteome data, some  were identified exclusively in the proteome dataset .
we observed that leukocyte related biological functions, such as activation of mononuclear leukocyte, mobilization of ca2+ and activation of lymphocyte were most strongly associated with the leukemia cell lines .
enrichment analysis suggested that the most significantly enriched pathways are, leukemia extravasation signaling pathway , which is responsible for leukocyte migration and related to metastasis in leukemia cell lines  <cit> , t cell receptor signaling  and icos-icosl signaling in t helper cells .
to further demonstrate the advantage of combining multiple layers of information in pathway analysis, we performed identical analysis only based on transcriptome markers from all of the four microarray studies.
although leukocyte extravasation signaling was still the most enriched pathway, it did not reach the same level of significance .
in addition, pathways that are not strongly associated with leukemia were also significantly enriched .
several pathways that are associated with leukemia and were detected in the combined analysis were absent, including nf-kb pathway and pi3k signaling in b lymphocytes .
we repeated this analysis on the set of mcia discovered features associated with melanoma .
the selected features comprised of proteins and genes that are highly expressed in melanoma cell lines, such as tyr, tyrp <dig> and bcl2a <dig>  these were significantly enriched in the biological functions or pathways associated with eumelanin biosynthesis and disorder of pigmentation including the melanocyte development and pigmentation signaling pathway .
melanocytic development and pigmentation is regulated in large part by the bhlh-lz microphthalmia-associated transcription factor  and mitf activity is controlled by at least two pathways: msh and kit signaling.
bcl2a <dig> is transcriptionally activated by mitf and serves as an anti-apoptosis factor  <cit> .
interestingly, the upstream regulator of mitf, lef <dig>  was also consistently identified on the same direction in all transcriptome datasets .
it is of note that, although all five datasets contributed to the coverage of this pathway, mitf was solely detected in the affymetrix data.
mcia can increase coverage and, the power of pathway  analyses as it does not require mapping or pre-filtering of features to the subset common to all datasets.
mcia allows easy integration of multiple omics levels to identify classes that are relevant in the given biological context.
comparison of mcia and regularized generalized canonical correlation analysis
in generalized canonical correlation analysis  several sets of variables are analyzed simultaneously.
several generalizations of cca have been described.
these employ different methods, including sum of correlations , sum of squared correlations  and sum of absolute value correlations   <cit> .
recently tenenhaus and coworkers introduced regularized generalized canonical correlation analysis  to generalize rcca to multi-block data analysis of data where the number of variables exceed the number of cases  <cit> .
we compare mcia to several rgcca methods that are defined by different shrinkage parameters and optimization criteria .
first, we compared three different optimization criteria in rgcca, namely sumcor, sabscor, ssqcor with mcia.
as depicted in additional file 1: figure s <dig>  the sumcor method and mcia algorithm consistently return similar results with positively correlated axes .
also the identified components from the sabscor and ssqcor methods are always highly correlated to the mcia results, but it is important to note that the correlation could be either positive or negative.
this is inconvenient for the comparison and integration of multiple omics datasets, as the components from one dataset might be inverted in another dataset.
by tuning the shrinkage parameter т, which can range from  <dig> to  <dig>  rgcca balances optimizing the intra-table and inter-table covariance.
additional file 1: figure s <dig> and s <dig> show that the identified components are nearly identical across datasets for т =  <dig>  the smaller the shrinkage parameter т, the higher is the correlation between neighboring components from different datasets.
but the variance of each individual dataset is less well explained by the components.
in contrast, the results of rgcca with a shrinkage parameter of т =  <dig> are very similar to mcia results.
in this case, rgcca gives priority to finding a component that explains its own block well  <cit> .
similarly, mcia maximizes the variance within each table and the covariance of components of each table with a consensus reference structure through a synthetic analysis.
it is important to note that in omics data analyses, the number of features is generally much larger than the number of observations.
therefore, a low т should be avoided as it results in overfitting of the data and apparently perfect correlations, which rarely represent meaningful information.
integrated analysis of microarray and rna-sequencing ovarian cancer datasets
in the ovarian cancer datasets, mcia was applied to several microarray and rna-seq gene expression datasets; agilent, affymetrix, rnaseq, rnaseqv <dig> which contained  <dig> ,  <dig> ,  <dig> , and  <dig>  genes respectively.
in the mcia space, the first pc  accounted for  <dig> % of the total variance and the second pc  accounted for  <dig> % of variance .
in comparison to microarray data, rna sequencing data typically contains many missing values.
these are generated when multiple experiments are combined.
we excluded genes  with high number of missing values.
after filtering genes with more than  <dig> missing values in rna-seq data, the four datasets contributed similarly to the total variance .
among the two rna-seq datasets, rnaseq consistently tended to be more variant than rnaseqv <dig> on pc1- <dig> .
rnaseq and rnaseqv <dig> were generated from the same illumina rna-sequencing data but using two different pre-processing approaches.
mcia results indicated that normalization and quantification with the rpkm method  outperforms mapsplice and rsem .
the informativeness or variance in rna sequencing data tended to be sensitive to pre-processing and filtering algorithms which is expected given that methods for processing these data are still emerging.
in addition, affymetrix profiles were generally more variant than agilent as indicated by greater pseudo-eigenvalues on pc1- <dig>  when the microarray and rnaseq data were compared, we detected several outlier genes that were highly variant on pc <dig> and pc <dig> on rnaseq but absent on the microarray platforms.
these include cdhr <dig> and hesrg which are highly expressed by the differentiated subtype   <cit> .
mcia identified ovarian subtypes
we applied mcia to compare the consistency and discrepancy in gene expression profiles of ovarian cancer tumors obtained by rna-sequencing and affymetrix and agilent microarray technologies .
the results revealed high overall similarity in structure between the four datasets and three platforms.
recent microarray gene expression profiling studies have reported four subtypes of ovarian cancer   <cit> .
these hgs-ovca subtypes can be clearly distinguished along the first two mcia axes .
the first axis generally separated samples with immunoreactive versus proliferative characteristics.
whereas the second axis distinguished tumors with a mesenchymal subtype which show a short survival time  <cit>  from the differentiated ovarian cancer samples.
consistent with other studies, mcia identified large overlap between the four subtypes, indicating that most samples exhibited features of multiple subtype signatures  <cit> .
in order to find whether this classification correlates with clinical factors, we compared the first two pcs with clinical records provided from the tcga data portal and the verhaak study  <cit> .
this comparison revealed that age at diagnosis was significantly negatively correlated with pc <dig> and positively correlated pc <dig> , suggesting that differentiated and immunoreactive patients tend to present at younger age.
the percentage of stromal cells is positively correlated with pc <dig> , which is in consensus with the mesenchymal subtype having greater percentage of stromal cells  <cit> .
other clinical factors, such as somatic mutation, drug treatment and tumor stages did not significantly correlate with either axis.
mcia suggests robust subtype biomarkers
both microarray and rna sequencing data resulted in a similar ordination of tumor samples in the mcia space.
in order to identify which genes contribute significantly to the divergence of samples, we examined the gene expression variables superimposed onto the same space .
the top  <dig> genes from the end of each axis were selected.
the full list of selected genes and their enriched pathways are provided in additional file 5: table s <dig>  each dataset contained different genes.
approximately 47% of genes  were measured on all four datasets .
among  <dig> genes selected as the top  <dig> genes from each datasets on pc <dig> and pc <dig> only  <dig> genes were in at least three platforms and  <dig>  were present in all datasets.
several of these “robust” markers, have been previously implicated in ovarian cancer  <cit> .
many t-cell activation and trafficking genes, such as cxcl <dig>  cd <dig> and cd3d were projected onto the positive end of the first axes, which represented the immunoreactive subtype tumors.
mcia revealed new markers that might be associated with the immune system, including sh2d1a, rhoh, saa <dig>  saa <dig> and gnly.
this is further corroborated by numerous go terms significantly associated with genes on this end of the axis   <cit> .
for instance, significantly enriched gene sets include glycoprotein, chemotaxis, defense and immune response .
the genes at the opposite end of the mcia axes included transcriptional factors sox <dig>  hmga <dig>  along with several cell cycle promoters, such as bex <dig>  mapk <dig> as well as nerve system development regulators , which characterize the proliferative subtype.
genes that are expressed on the positive end of axis  <dig>  such as postn, cxcl <dig> and hoxa <dig>  define the mesenchymal cluster.
other potential mesenchymal subtype markers for ovarian cancer include aspn, homeobox alpha genes as well as collagens.
aspn is a critical regulator of tgf-beta pathway that induces the epithelial mesenchymal transition.
gene set analysis revealed that mesenchymal genes are enriched in go terms including cell adhesion, skeletal system development, collagen and ecm receptor interaction pathway .
the robust markers at the differentiated end include oviductal glycoprotein  <dig> , spdef, kiaa <dig>  gjb <dig> and alppl <dig>  some of which have already been reported as ovarian biomarkers.
for instance, ovgp <dig> has been suggested as a possible serum marker for the detection of low grade ovarian cancer  <cit> .
although the tcga dataset is all high grade serous ovarian cancer, in our analysis, it was highly expressed in differentiated subtype.
human spdef protein plays a significant role in tumorigenesis in multiple cancers, including ovarian cancer and has been reported to suppress prostate tumor metastasis.
a recent study on prostate cancer demonstrated that spdef suppresses cancer metastasis through down-regulation of matrix metalloproteinase  <dig> and  <dig> , which are required for the invasive phenotype of cells  <cit> .
our analysis implied that spdef and matrix metalloproteinase plays a similar role in the development of ovarian cancer.
in addition, it has been shown that, in a mouse model, postn down-regulates alpp mrna  <cit> .
postn and alppl <dig> were projected onto the diametral ends of axes  <dig>  which implies that the same mechanism of regulation exists in ovarian cancer and can be exploited to distinguish subtypes.
interestingly, the david gene set analysis of markers for the differentiated phenotype did not reveal as strong gene set enrichments as described for the other subtypes  indicating that this subtype exhibits considerably higher degree of heterogeneity.
CONCLUSIONS
in the present study, we described multiple co-inertia analysis , an exploratory data analysis method that can identify co-relationships between multiple high dimensional datasets.
mcia projects multiple sets of features onto the same dimensional space and provides a simple graphical representation for the efficient identification of concordance between datasets.
the sets of features may have none or few features in common.
by transforming multiple sources of data onto the same scale, the most variant features are transformed onto the same scale.
this allows one to extract and easily combine sets of omic features  for greater power in subsequent pathway analysis.
mcia provides a consensus reference structure of datasets, revealing similar trends among multiple tables.
compared to rgcca, we found that mcia is most similar to the sumcor version of rgcca with т =  <dig> in practice.
our integrative analysis of nci- <dig> cell line panel indicated that, although both transcriptome and proteome cell lines were clustered according to their lineage, they provides complementary information.
we demonstrated that integrated analysis of gene and protein expression data increases the power of pathway analysis and yields more information than an analysis of gene expression alone.
mcia highly ranked the leukemia extravasation signaling pathway.
this pathway was overrepresented with features that were predominantly from the proteomics data and were enriched in biological functions of “activation of mononuclear leukocyte and lymphocyte”.
mcia of high grade serous ovarian cancer revealed four previously described subtypes of ovarian cancer and provided novel subtype markers.
an advantage of mcia is that it couples multiple set of features measured on the same set of samples.
since it does not rely on feature annotation, it is not limited by the immaturity of annotations.
there is no prerequisite to filter or map features  to a common set thereby considerably increasing genome coverage.
in a study that compares cia with other sparse multiple table analysis methods , lecao et al.
suggested that cia may result in redundancy when it is used for feature selection since it does not include a built-in procedure for variable selection  <cit> .
similarly, mcia does not impose any sparsity in the result, so mcia selects much more features than methods introducing the lasso penalty, such as sgcca  <cit>  or pcca  <cit> .
hence, the interpretation of mcia selected features would have to be coupled with other methods, such as enrichment analysis, in order to reveal functional insights.
we also note that the mcia algorithm finds solutions in a sequential manner and each order of components requires a singular value decomposition  for a large dataset.
the computationally intensity of the algorithm increases with sample size as more components are retained.
for instance, the cpu time of analysis of the nci- <dig> data with  <dig> kept principal components was around  <dig> seconds on intel xeon  <dig> mhz using one thread of a linux server.
in conclusion, we believe mcia is a useful method for integration of multiple omics datasets where the same tissue or cell lines have been assayed multiple times.
mcia is available to the community via an r-bioconductor  package which includes documentation and a vignette.
availability of supporting data
the microarray data of nci- <dig> cell lines are available through cellminer .
the nci- <dig> proteomic data can be downloaded from http://wzw.tum.de/proteomics/nci60/ as well as from https://www.proteomicsdb.org.
the ovarian cancer data are available through the tcga download portal .
abbreviations
cca: canonical correlation analysis; cca-en: canonical correlation analysis with elastic net penalty; cia: co-inertia analysis; cns: central nervous system; coa: correspondence analysis; cpca: consensus principal component analysis; emt: epithelial-mesenchymal transition; encode: the encyclopedia of dna elements; gcca: generalized canonical correlation analysis; gcrma: gc robust multichip averaging; gelc-ms/ms: in-gel digestion and liquid chromatography tandem mass spectrometry; hgs-ovca: high grade serous ovarian cancer; ica: independent component analysis; ipa: ingenuity pathways analysis; ipi: international protein index; mcia: multiple co-inertia analysis; ms: mass spectrometry; nci: the national cancer institute; nsc: non-symmetric correspondence analysis; nsclc: non-small-cell lung carcinoma; page: polyacrylamide gel electrophoresis; pc: principal component; pca: principal component analysis; pcca: penalized canonical correlation analysis; rgcca: regularized generalized canonical correlation analysis; rma: robust multichip averaging; rnaseq: rna sequencing; rpkm: reads per kilo base per million; rsem: rna-seq by expectation maximization; svd: singular value decomposition; sgcca: sparse generalized canonical correlation analysis; spls: sparse partial least square; tcga: the cancer genome atlas.
competing interests
the authors declare that they have no competing interests.
authors’ contribution
cm carried out the analysis as a phd student in the group of bk and wrote the manuscript.
bk provided input regarding the interpretation of the results.
ac and amg made numerous important intellectual contributions, provided input for both the design of the study and drafting of the manuscript.
amg designed and supervised the study and wrote the manuscript.
all authors read and approved the manuscript.
supplementary material
additional file 1
supplementary information.
figures s1–s <dig>
click here for file
additional file 2: table s1
full list of biological markers highly weighted on each mcia axis of the nci- <dig> data.
click here for file
additional file 3: table s2
full list of leukemia and melanoma markers  and the ipa functional analysis.
click here for file
additional file 4: table s3
pathway analysis of leukemia markers.
click here for file
additional file 5: table s4
functional analysis of the ovarian subtype markers.
click here for file
acknowledgements
the authors gratefully acknowledge hannes hahne and mathias wilhelm of tum and curtis huttenhower of harvard school of public health for critical reading of the manuscript.
funding for this work was provided by dana-farber cancer institute women’s cancers program, the claudia adams barr foundation and the national cancer institute at the national institutes of health .

 there has been rapid progress in high-throughput technologies and platforms to assay global mrna, mirna, methylation, proteins, and metabolite profiles of cells are readily available. advances in rna-sequencing and mass spectrometry based proteomics have dramatically improved coverage and quality of genomic, transcriptomic and proteomic profiling. increasing number of studies including the cancer genome atlas and encyclopedia of dna elements projects systematically profile large number of biological samples resulting in multiple levels of quantitative information. recent advances of ms based proteomics provide a complementary approach to genomics and transcriptomic technologies and systematic analyses can now be carried out to identify and quantify the majority of proteins expressed in human cells. these data yield unprecedented views of molecular building blocks and the machinery of cells. interpreting these large-scale datasets to derive information about a biological system represents a considerable challenge often faced by investigators. multiple omics data analysis can be broadly defined by some common questions, which are dependent on the data collected; multiple datasets measuring the same biological molecules or multiple datasets each measuring different biological molecules. in the first case, given multiple transcriptomics data from different microarray or rna-sequencing studies, the aim may be to discover which platform is the most informative with highest quality data, identify robust biomarkers across datasets or highlight platform specific discrepancies in measurements. in the second case, given multiple different data such as transcripts, proteins and metabolites, the objective may be to integrate and concatenate information to increase the breath and coverage of available data in a biological network. in this case, specific platform discrepancies are less important and performance of data integration is more likely to be assessed using system biology or pathway approaches. nevertheless, both analyses face common challenges associated with integrating data from disparate technologies. several meta-analysis studies map identifiers from each platform to a common set of identifiers to generate a single concatenated matrix for subsequent analysis. however, this data simplification overlooks several fundamental platform and biological biases. platforms are not universal and measure different molecules. filtering genes to their intersection may considerably reduce data coverage. in addition, the many-to-many mapping of gene identifiers from multiple platforms complicates direct comparison of molecules across multiple levels. moreover, because correlations between different platforms are probably lower than expected, it may not provide gains in data quality or study power. such filtering may also introduce bias because platform discrepancies could reflect biological variation. for instance, poor correlation between a transcript and its translated protein may result from biological processes such as microrna post-transcriptional repression. similarly, correlations between proteins and metabolites of pathways can diverge if proteins are expressed in an inactive form, in which case its abundance may not represent activity. ordination methods, such as principal component analysis, independent component analysis and correspondence analysis, are exploratory data analysis approaches that have been applied to analyze omics data including transcriptome and proteome studies. graphical representation of measurements and variables on a lower dimensional space facilitates interpretation of global variance structure and identification of the most informative features across datasets. these methods permit visualization of data that have considerable levels of noise and data where the number of variables exceeds the number of measurements, which is typical in omics studies. however, these approaches do not solve the problem of comparing many datasets simultaneously. studies have extended these approaches to couple two datasets together. one such approach is co-inertia analysis. cia was originally applied to study ecological and environmental tables, where it was employed to link environmental variables with species characteristics. culhane and colleagues introduced cia in genomics, when they compared data from two microarray platforms. an advantage of this method is that it does not require the mapping or filtering of genes to a common set. the relationship between co-inertia analysis and related methods including procrustes analysis, canonical correlation analysis with elastic net penalization and sparse partial least squares have been described previously. cia and spls both maximize the covariance between eigenvectors and are efficient in determining main individual effects in paired dataset analysis. by contrast cca-en maximizes the correlation between eigenvectors and tends to discover effects present in both datasets, but may omit to discover strong individual effects. variables selected by cca-en and spls are highly similar but cia selected marginally different marker genes that may have some redundancy. a noteworthy advantage of cia is that it can be coupled with several dimension reduction approaches, including pca or correspondence analysis, such that it can accommodate both discrete count data and continuous data. these approaches are performed on each dataset separately and can be integrated using cia. however, all above methods including cia are limited to the analysis of two datasets, limiting their application in modern multi-omics studies. several approaches have been proposed for integrating more than two datasets, such as consensus pca, regularized generalized canonical correlation analysis, sparse generalized canonical correlation analysis and penalized canonical correlation analysis. sgcca and pcca originally focus on the feature selection from multiple datasets, but also can be used for multiple table integration problem. here, we describe another method, multiple co-inertia analysis, for the analysis of more than two omics datasets, extending its application in the field of environmental science and, recently, phylogenetics. mcia is related to consensus pca which both maximize the square covariance between eigenvectors and are subject to similar constraints. cpca is less sensitive to multi-collinearity within each dataset than generalized canonical correlation analysis. we illustrate the application of mcia using two different examples, and show that integrated analysis is more insightful than analysis of the individual datasets. first, we demonstrated the power of mcia via applying it to the integration and comparison of multi-omics data independent of data annotation. we employed mcia to identify common relationships among multiple gene and protein expression data of the nci- cancer cell line panel of the national cancer institute. the integrated analysis revealed that cell lines are clustered according to anatomical tissue source and showed a significant degree of correlation between transcript and protein expression. second, we assessed the concordance in gene expression data obtained from microarray and next generation rna-sequencing of samples of high grade serous ovarian cancer. mcia integrated ovarian cancer gene expression data from different sources which captured distinct subsets of the transcriptome to reveal a set of biomarkers that were consistently highly ranked by all four platforms and were biologically relevant to ovarian cancer. to enable community access to mcia, we implemented the method into the r-bioconductor package as an easy-to-use tool for bioinformaticians and biologists. methods mathematical basis of mcia a typical omics dataset is a matrix where the number of features exceeds the number of measurements. mcia requires a set of tables where either features or measurements are matched and have equal weights. mcia is performed in a two-step process. first a one table ordination method, such as pca, coa or non-symmetric correspondence analysis is applied on each dataset separately, which transforms data into comparable lower dimensional spaces. in our analysis, given an omics data table m =  with 1 ≤ i ≤ n and 1 ≤ j ≤ q, where m is a matrix, i indicates row index and j for column index. we denote the row and column sums of m as mi+ and m+j respectively, and m++ as the grand total. the relative contribution or weight of row i to the total variation in the data set is denoted ri and calculated as ri = mi+/m++ while the relative contribution of column j is denoted as cj = m+j/m++. similarly, the contribution of each individual element of m to the total variation pij can be calculated as pij = mij/m++. we then derive a new matrix x with the values defined above as xij=pijri−cj where xij is the centered row profile, i.e. the relative abundance of selected variable to the measurement’s weight. the second step in mcia is a generalization of cia. it solves the problem of simultaneous analysis of a set of statistical triplet where k = …, k,…, k and xk is a set of transformed matrices. qk is a qk × qk matrix with rij in diagonal elements, indicating the hyperspace of features metrics. d is an n × n matrix which is an identity matrix indicating equal weight across all columns in all tables. mcia maximizes the sum of the squared covariance between scores of each table with synthetic axes ν, that is: fu …,uk,…,uk,v=∑k=1kwkcov2xkqkuk,v where cov stands for the square of covariance of quantities inside parenthesis and ωk is the weight of each table. the v represents the reference structure or synthetic center and uk are auxiliary axes. the score of each individual table would then be xkqkuk. in contrast with other ordination methods, mcia finds solutions sequentially. multiple matrices xk can be weighted and concatenated to a single matrix x =. similarly, a single feature metric q could be concatenated as q =.the first order solutions of u to uk and v are given by the first principal component of the following eigen-system: wxqxtdv=λv then the normalized auxiliary axis uk are uk1=xktdv1/||xktdv1||qkk= …,k where ||•|| is the norm in the qk metric. the subsequent solutions are found with residual matrices from the calculation of the first order solution with the constraint that the remaining order axes are orthogonal with the previous sets, namely: vjtdvs=0andukjtqkuks=01≤j<s the residual matrices used by second order solution is deflated as x1order2=x1‒x1pk where the projection matrix pk is pk1=uk1uk1qkuk1t‒1uk1qk the superscript t and − stand for matrix transposition and matrix inversion respectively. therefore, the formula removes the dimension that is spanned by vector uk to get a residual matrix, which is passed to the svd to find the second order solution. these steps are repeated until the desired number of axes is generated. as a result, mcia provides a simultaneous ordination of columns and rows of multiple tables within the same hyperspace, with features or measurements sharing similar trends will be closely projected. the detailed description of mcia and the proof that these axes are maximally co-variant are given in chessel and hanafi. datasets we analyzed publicly available sets of data from two studies: transcriptomic and proteomic datasets of the nci- cancer cell line panel, the latter one generated in our group, and an ovarian cancer dataset generated as part of the tcga project. in each study, there are multiple datasets measuring molecules from the same samples. nci- data the nci- panel is a collection of cancer cell lines of leukemia, lymphomas, melanomas and carcinomas of ovarian, renal, breast, prostate, colon, lung and central nervous system origin. the nci- transcriptome data were downloaded from cellminer and were obtained on four different platforms; affymetrix hg-u plus, hg-u hg-u and agilent ge 4x44k. affymetrix data were normalized using gc robust multichip averaging gcrma; and agilent data were log transformed as obtained from the cellminer. although data filtering is not required to perform mcia, to facilitate data interpretation, microarray data were filtered to exclude probes that do not map to an official hugo gene symbol. the probe with highest average value was retained when multiple probes mapped to the same gene. filtering produced datasets of;; and genes on agilent, hg-u hg-u and hg-u plus platforms respectively. the lung cancer cell line nci-h was excluded since its expression profile was not available on the hg-u platform. a venn diagram representing the overlapping genes in the processed data for each platform is provided in additional file 1: figure s the proteome profiles of cell lines were produced from a conventional gelc-ms/ms approach and label-free quantification, as described in. the international protein index identifiers were mapped to official gene symbol to facilitate subsequent pathway interpretation. data were log transformed and no filtering or additional normalization were performed. this dataset represents protein expressions across cell line in nci- panel. ovarian cancer datasets gene expression of tumors from ovarian cancer patients were profiled using two microarray platforms and rna-sequencing on illumina hiseq platform. data were downloaded from the nci-tcga data portal 07/08/2013;. patient samples that were present in all four datasets were included in the analysis. the agilent and affymetrix data were normalized and summarized by lowess and robust multichip averaging, respectively. the transcript expression levels of the illumina rna-sequencing data were determined using two different pre-processing pipelines denoted as rnaseq and rnaseqv respectively. normalization and quantification of rnaseq followed the rpkm method whereas the alignment and gene expression quantification in rnaseqv were obtained by mapsplice and rsem. in rnaseq and rnaseqv2; and genes were detected. these data were filtered to exclude genes with more than missing values. only genes mapped to an official gene symbol were retained. for the features mapped to the same gene symbol, the one with the largest average expression value was kept. remaining missing values were replaced with a positive value far smaller than the lowest expression value in each dataset and then, the expression values were log transformed. after filtering, the agilent, affymetrix, rnaseq and rnaseqv datasets contained;; and gene expression measurements respectively. the venn diagram representing the overlap of genes in these datasets is shown in additional file 1: figure s  integrated analysis of the nci- cell line transcriptome and proteome the nci- panel, a collection of cancer cell lines derived from nine different tissues has been extensively used in in vitro high-throughput drug screen assays. they have been molecularly profiled using comparative genomic hybridization array, karyotype analysis, dna mutational analysis, transcripts expression array, microarrays for microrna expression and protein expression. mcia was applied as an exploratory analysis of four transcriptomic studies and one proteomic study of the cell lines. figure 1a shows the projection of cell lines onto the first two principal components of mcia. similar to the visualization employed in cia, the datasets are transformed into the same projection. the coordinates of the five measurements for each cell line are connected by lines. the length of which indicates the divergence between the mrnas and protein expression levels for a particular cell line. the mcia plot of the first two principal components shows similar trends in transcriptome and proteome profiles, indicating that the most variant sources of biological information were similar. cell lines originating from the same or closely related anatomical source of tissue were projected close to each other and converged into groups. the colon, leukemia, melanoma, cns, renal and ovarian cell lines segregated largely according to their tissue of origin. seven out of eight melanoma lines clustered together, and the remaining one, lox-imvi, has been reported to lack melanin production. these results are consistent with independently performed hierarchical clustering analysis. there was greater divergence in the cell lines from tumors with more intrinsic molecular heterogeneity. the transcriptome and proteome profiles of the individual breast and lung cell lines were projected close in space demonstrating that the expression profiles shared a high degree of consensus. the tight projection of multiple data types from diverse technology platforms provides evidence that the observed spread of cell lines reflected the biological variance, as opposed to inter-study technical or other stochastic variance. for instance, we observed that the estrogen receptor positive breast cancer cell line mcf displays an epithelial phenotype and clustered to colon cancer lines. in contrast, the cell line negative for the estrogen receptor, hs578t, clustered with the stromal/mesenchymal cluster of glioblastoma and renal tumor cell lines. this suggests that hs578t exhibits more invasive mesenchymal features compared to t47d and mcf overall co-structure comparison using mcia each pc has an associated eigenvalue which represents the amount of variability contained in that pc. the first three pcs of the mcia accounted for %, % and % of variance respectively. the observation that the first two pcs capture less than a third of the structure in the datasets reflects the complexity inherent in cell lines of tumors from nine different organs. in order to identify the contribution of each dataset to the total variance, that is, to what extent each dataset deviates or agrees with what the majority of datasets support, we extracted the mcia pseudo-eigenvalues. figure 1b shows the pseudo-eigenvalues associated with the first two principal components of each dataset. comparison within microarray data revealed that affymetrix hgu plus accounts for the highest variance on axis and possibly because this platform contains informative features, or features that are poorly detected or absent on other platforms. we observed that the similarity within transcriptome datasets is greater than the similarity between transcriptome and proteome data, which is consistent with the results shown by figure 1a. one of the most attractive features of mcia is that it can be used to highlight lack or presence of co-structure between datasets, thus it allows selection of the strongest features from each dataset for subsequent analysis. for instance, we observed in particular, large variation between the protein and transcript expression patterns of two cell lines, melanoma skmel and ovarian igrov the proteome coordinates of skmel were close to the origin and far from the transcriptomic data that was projected on the negative end of pc with the other melanoma cell line data. the poor information content in proteome data of the skmel cell line could reflect the lack of expression of melanin related genes on protein level. similarly, the incongruence of the proteome and transcriptome data of the ovarian cell line igrov may be due to expression of less epithelial markers that projected on the positive direction of axis to characterize the overall correlation between each pair of high dimensional data we calculated the pair-wise rv coefficient, a multivariate generalization of the squared pearson correlation coefficient. for each pair of datasets, the rv-coefficient is calculated as the total co-inertia divided by the square root of the product of the squared total inertia from the individual analysis. as the co-structure between two datasets increases, the rv score move towards to a zero rv score indicates no similarity. the overall similarity in structure between microarray data was higher than the similarity between microarray and proteomics data; average rv coefficient > and respectively. when mcia was performed on the same transcriptome data and the subset of proteome data that were quantified in all cell lines, the filtered proteome data had a higher consensus to the co-structure and increased pseudo-eigenvalues. mcia axes describe biological properties in contrast to traditional clustering methods, mcia projects the original data onto a lower dimensional space, maximizing the covariance of each dataset with respect to the reference structure. in mcia plots, a gene that is highly expressed in a certain cell line will be projected in the direction of this cell line and the greater the distance from the origin, the stronger the association. in order to identify biomarkers that are highly associated with cancer cell lines of different origins, we examined the feature space of mrnas and proteins that were projected in the same direction and space. the first axis, which explains the largest variance, separated cells with epithelial or mesenchymal characteristics, suggesting that epithelial-mesenchymal transition is an essential mechanism defining different classes of cancers. emt has been shown to play an important role in epithelial cell malignancy and metastasis. epithelial markers, including slc27a cdh spint s100p and epcam had high weights on the positive side of pc. at the opposite end, mesenchymal and collagen markers, including gja which is involved in epicardial to mesenchymal cell transition, and tgfβ were observed. the second axis, pc clearly separated melanoma and leukemia from other epithelial cancer types. the strongest determinant of the vertical axis is a set of melanoma-related genes, namely melanoma-associated antigens, melanogenic enzyme as well as tyr, dct, tyrp malana, s100b and bcl2a the top genes with greatest weights on pc and pc were selected from each dataset and the full list of markers is provided in additional file 2: table s among selected genes, were measured in three or more datasets. mcia enables the study of the union of features from all studies. among the nci- datasets, less than 12% of the total genes studied were measured in all five datasets. by observing highly ranked genes across studies, one can identify robust markers that could be subject to further analysis. integration of proteomics and transcriptomics complements the biological information to further evaluate the biological significance of the features selected by mcia, we employed ingenuity pathway analysis to discover significant canonical pathways which discriminate different cell lines. in mcia plots, samples and features are projected onto the same space. the genes with strongest association to a cell line are those projected in the same direction and have the highest weights. as features have been transformed on the same scale, the union of features from each individual dataset can be easily extracted and concatenated to provide greater coverage in pathway analysis. features strongly associated to each tissue type from both transcriptome and proteome datasets can be concatenated and mapped to signaling pathways. there is no requirement to extract equal numbers of features from each data type. for example we observed that features strongly associated with leukemia related features tended to be enriched in the proteins. the most extreme features associated with the leukemia cell lines were selected from all platforms using their coordinates and were subjected to the functional and pathway analysis. the full list of features, the coordinate feature selection criteria and their functional and pathway analysis are provided in additional file 3: tables s and additional file 4: table s complementary information can be obtained by integrating data from different platforms and data types which increases the genome coverage and power of subsequent pathway analysis. while numerous genes were over-expressed in both the transcriptome and proteome data, some were identified exclusively in the proteome dataset. we observed that leukocyte related biological functions, such as activation of mononuclear leukocyte, mobilization of ca2+ and activation of lymphocyte were most strongly associated with the leukemia cell lines. enrichment analysis suggested that the most significantly enriched pathways are, leukemia extravasation signaling pathway, which is responsible for leukocyte migration and related to metastasis in leukemia cell lines, t cell receptor signaling and icos-icosl signaling in t helper cells. to further demonstrate the advantage of combining multiple layers of information in pathway analysis, we performed identical analysis only based on transcriptome markers from all of the four microarray studies. although leukocyte extravasation signaling was still the most enriched pathway, it did not reach the same level of significance. in addition, pathways that are not strongly associated with leukemia were also significantly enriched. several pathways that are associated with leukemia and were detected in the combined analysis were absent, including nf-kb pathway and pi3k signaling in b lymphocytes. we repeated this analysis on the set of mcia discovered features associated with melanoma. the selected features comprised of proteins and genes that are highly expressed in melanoma cell lines, such as tyr, tyrp and bcl2a these were significantly enriched in the biological functions or pathways associated with eumelanin biosynthesis and disorder of pigmentation including the melanocyte development and pigmentation signaling pathway. melanocytic development and pigmentation is regulated in large part by the bhlh-lz microphthalmia-associated transcription factor and mitf activity is controlled by at least two pathways: msh and kit signaling. bcl2a is transcriptionally activated by mitf and serves as an anti-apoptosis factor. interestingly, the upstream regulator of mitf, lef was also consistently identified on the same direction in all transcriptome datasets. it is of note that, although all five datasets contributed to the coverage of this pathway, mitf was solely detected in the affymetrix data. mcia can increase coverage and, the power of pathway analyses as it does not require mapping or pre-filtering of features to the subset common to all datasets. mcia allows easy integration of multiple omics levels to identify classes that are relevant in the given biological context. comparison of mcia and regularized generalized canonical correlation analysis in generalized canonical correlation analysis several sets of variables are analyzed simultaneously. several generalizations of cca have been described. these employ different methods, including sum of correlations, sum of squared correlations and sum of absolute value correlations. recently tenenhaus and coworkers introduced regularized generalized canonical correlation analysis to generalize rcca to multi-block data analysis of data where the number of variables exceed the number of cases. we compare mcia to several rgcca methods that are defined by different shrinkage parameters and optimization criteria. first, we compared three different optimization criteria in rgcca, namely sumcor, sabscor, ssqcor with mcia. as depicted in additional file 1: figure s the sumcor method and mcia algorithm consistently return similar results with positively correlated axes. also the identified components from the sabscor and ssqcor methods are always highly correlated to the mcia results, but it is important to note that the correlation could be either positive or negative. this is inconvenient for the comparison and integration of multiple omics datasets, as the components from one dataset might be inverted in another dataset. by tuning the shrinkage parameter т, which can range from to rgcca balances optimizing the intra-table and inter-table covariance. additional file 1: figure s and s show that the identified components are nearly identical across datasets for т = the smaller the shrinkage parameter т, the higher is the correlation between neighboring components from different datasets. but the variance of each individual dataset is less well explained by the components. in contrast, the results of rgcca with a shrinkage parameter of т = are very similar to mcia results. in this case, rgcca gives priority to finding a component that explains its own block well. similarly, mcia maximizes the variance within each table and the covariance of components of each table with a consensus reference structure through a synthetic analysis. it is important to note that in omics data analyses, the number of features is generally much larger than the number of observations. therefore, a low т should be avoided as it results in overfitting of the data and apparently perfect correlations, which rarely represent meaningful information. integrated analysis of microarray and rna-sequencing ovarian cancer datasets in the ovarian cancer datasets, mcia was applied to several microarray and rna-seq gene expression datasets; agilent, affymetrix, rnaseq, rnaseqv which contained,,, and genes respectively. in the mcia space, the first pc accounted for % of the total variance and the second pc accounted for % of variance. in comparison to microarray data, rna sequencing data typically contains many missing values. these are generated when multiple experiments are combined. we excluded genes with high number of missing values. after filtering genes with more than missing values in rna-seq data, the four datasets contributed similarly to the total variance. among the two rna-seq datasets, rnaseq consistently tended to be more variant than rnaseqv on pc1-. rnaseq and rnaseqv were generated from the same illumina rna-sequencing data but using two different pre-processing approaches. mcia results indicated that normalization and quantification with the rpkm method outperforms mapsplice and rsem. the informativeness or variance in rna sequencing data tended to be sensitive to pre-processing and filtering algorithms which is expected given that methods for processing these data are still emerging. in addition, affymetrix profiles were generally more variant than agilent as indicated by greater pseudo-eigenvalues on pc1- when the microarray and rnaseq data were compared, we detected several outlier genes that were highly variant on pc and pc on rnaseq but absent on the microarray platforms. these include cdhr and hesrg which are highly expressed by the differentiated subtype. mcia identified ovarian subtypes we applied mcia to compare the consistency and discrepancy in gene expression profiles of ovarian cancer tumors obtained by rna-sequencing and affymetrix and agilent microarray technologies. the results revealed high overall similarity in structure between the four datasets and three platforms. recent microarray gene expression profiling studies have reported four subtypes of ovarian cancer. these hgs-ovca subtypes can be clearly distinguished along the first two mcia axes. the first axis generally separated samples with immunoreactive versus proliferative characteristics. whereas the second axis distinguished tumors with a mesenchymal subtype which show a short survival time from the differentiated ovarian cancer samples. consistent with other studies, mcia identified large overlap between the four subtypes, indicating that most samples exhibited features of multiple subtype signatures. in order to find whether this classification correlates with clinical factors, we compared the first two pcs with clinical records provided from the tcga data portal and the verhaak study. this comparison revealed that age at diagnosis was significantly negatively correlated with pc and positively correlated pc, suggesting that differentiated and immunoreactive patients tend to present at younger age. the percentage of stromal cells is positively correlated with pc, which is in consensus with the mesenchymal subtype having greater percentage of stromal cells. other clinical factors, such as somatic mutation, drug treatment and tumor stages did not significantly correlate with either axis. mcia suggests robust subtype biomarkers both microarray and rna sequencing data resulted in a similar ordination of tumor samples in the mcia space. in order to identify which genes contribute significantly to the divergence of samples, we examined the gene expression variables superimposed onto the same space. the top genes from the end of each axis were selected. the full list of selected genes and their enriched pathways are provided in additional file 5: table s each dataset contained different genes. approximately 47% of genes were measured on all four datasets. among genes selected as the top genes from each datasets on pc and pc only genes were in at least three platforms and were present in all datasets. several of these “robust” markers, have been previously implicated in ovarian cancer. many t-cell activation and trafficking genes, such as cxcl cd and cd3d were projected onto the positive end of the first axes, which represented the immunoreactive subtype tumors. mcia revealed new markers that might be associated with the immune system, including sh2d1a, rhoh, saa saa and gnly. this is further corroborated by numerous go terms significantly associated with genes on this end of the axis. for instance, significantly enriched gene sets include glycoprotein, chemotaxis, defense and immune response. the genes at the opposite end of the mcia axes included transcriptional factors sox hmga along with several cell cycle promoters, such as bex mapk as well as nerve system development regulators, which characterize the proliferative subtype. genes that are expressed on the positive end of axis such as postn, cxcl and hoxa define the mesenchymal cluster. other potential mesenchymal subtype markers for ovarian cancer include aspn, homeobox alpha genes as well as collagens. aspn is a critical regulator of tgf-beta pathway that induces the epithelial mesenchymal transition. gene set analysis revealed that mesenchymal genes are enriched in go terms including cell adhesion, skeletal system development, collagen and ecm receptor interaction pathway. the robust markers at the differentiated end include oviductal glycoprotein, spdef, kiaa gjb and alppl some of which have already been reported as ovarian biomarkers. for instance, ovgp has been suggested as a possible serum marker for the detection of low grade ovarian cancer. although the tcga dataset is all high grade serous ovarian cancer, in our analysis, it was highly expressed in differentiated subtype. human spdef protein plays a significant role in tumorigenesis in multiple cancers, including ovarian cancer and has been reported to suppress prostate tumor metastasis. a recent study on prostate cancer demonstrated that spdef suppresses cancer metastasis through down-regulation of matrix metalloproteinase and, which are required for the invasive phenotype of cells. our analysis implied that spdef and matrix metalloproteinase plays a similar role in the development of ovarian cancer. in addition, it has been shown that, in a mouse model, postn down-regulates alpp mrna. postn and alppl were projected onto the diametral ends of axes which implies that the same mechanism of regulation exists in ovarian cancer and can be exploited to distinguish subtypes. interestingly, the david gene set analysis of markers for the differentiated phenotype did not reveal as strong gene set enrichments as described for the other subtypes indicating that this subtype exhibits considerably higher degree of heterogeneity.  in the present study, we described multiple co-inertia analysis, an exploratory data analysis method that can identify co-relationships between multiple high dimensional datasets. mcia projects multiple sets of features onto the same dimensional space and provides a simple graphical representation for the efficient identification of concordance between datasets. the sets of features may have none or few features in common. by transforming multiple sources of data onto the same scale, the most variant features are transformed onto the same scale. this allows one to extract and easily combine sets of omic features for greater power in subsequent pathway analysis. mcia provides a consensus reference structure of datasets, revealing similar trends among multiple tables. compared to rgcca, we found that mcia is most similar to the sumcor version of rgcca with т = in practice. our integrative analysis of nci- cell line panel indicated that, although both transcriptome and proteome cell lines were clustered according to their lineage, they provides complementary information. we demonstrated that integrated analysis of gene and protein expression data increases the power of pathway analysis and yields more information than an analysis of gene expression alone. mcia highly ranked the leukemia extravasation signaling pathway. this pathway was overrepresented with features that were predominantly from the proteomics data and were enriched in biological functions of “activation of mononuclear leukocyte and lymphocyte”. mcia of high grade serous ovarian cancer revealed four previously described subtypes of ovarian cancer and provided novel subtype markers. an advantage of mcia is that it couples multiple set of features measured on the same set of samples. since it does not rely on feature annotation, it is not limited by the immaturity of annotations. there is no prerequisite to filter or map features to a common set thereby considerably increasing genome coverage. in a study that compares cia with other sparse multiple table analysis methods, lecao et al. suggested that cia may result in redundancy when it is used for feature selection since it does not include a built-in procedure for variable selection. similarly, mcia does not impose any sparsity in the result, so mcia selects much more features than methods introducing the lasso penalty, such as sgcca or pcca. hence, the interpretation of mcia selected features would have to be coupled with other methods, such as enrichment analysis, in order to reveal functional insights. we also note that the mcia algorithm finds solutions in a sequential manner and each order of components requires a singular value decomposition for a large dataset. the computationally intensity of the algorithm increases with sample size as more components are retained. for instance, the cpu time of analysis of the nci- data with kept principal components was around seconds on intel xeon mhz using one thread of a linux server. in conclusion, we believe mcia is a useful method for integration of multiple omics datasets where the same tissue or cell lines have been assayed multiple times. mcia is available to the community via an r-bioconductor package which includes documentation and a vignette. availability of supporting data the microarray data of nci- cell lines are available through cellminer. the nci- proteomic data can be downloaded from http://wzw.tum.de/proteomics/nci60/ as well as from https://www.proteomicsdb.org. the ovarian cancer data are available through the tcga download portal. abbreviations cca: canonical correlation analysis; cca-en: canonical correlation analysis with elastic net penalty; cia: co-inertia analysis; cns: central nervous system; coa: correspondence analysis; cpca: consensus principal component analysis; emt: epithelial-mesenchymal transition; encode: the encyclopedia of dna elements; gcca: generalized canonical correlation analysis; gcrma: gc robust multichip averaging; gelc-ms/ms: in-gel digestion and liquid chromatography tandem mass spectrometry; hgs-ovca: high grade serous ovarian cancer; ica: independent component analysis; ipa: ingenuity pathways analysis; ipi: international protein index; mcia: multiple co-inertia analysis; ms: mass spectrometry; nci: the national cancer institute; nsc: non-symmetric correspondence analysis; nsclc: non-small-cell lung carcinoma; page: polyacrylamide gel electrophoresis; pc: principal component; pca: principal component analysis; pcca: penalized canonical correlation analysis; rgcca: regularized generalized canonical correlation analysis; rma: robust multichip averaging; rnaseq: rna sequencing; rpkm: reads per kilo base per million; rsem: rna-seq by expectation maximization; svd: singular value decomposition; sgcca: sparse generalized canonical correlation analysis; spls: sparse partial least square; tcga: the cancer genome atlas. competing interests the authors declare that they have no competing interests. authors’ contribution cm carried out the analysis as a phd student in the group of bk and wrote the manuscript. bk provided input regarding the interpretation of the results. ac and amg made numerous important intellectual contributions, provided input for both the design of the study and drafting of the manuscript. amg designed and supervised the study and wrote the manuscript. all authors read and approved the manuscript. supplementary material additional file 1 supplementary information. figures s1–s click here for file additional file 2: table s1 full list of biological markers highly weighted on each mcia axis of the nci- data. click here for file additional file 3: table s2 full list of leukemia and melanoma markers and the ipa functional analysis. click here for file additional file 4: table s3 pathway analysis of leukemia markers. click here for file additional file 5: table s4 functional analysis of the ovarian subtype markers. click here for file acknowledgements the authors gratefully acknowledge hannes hahne and mathias wilhelm of tum and curtis huttenhower of harvard school of public health for critical reading of the manuscript. funding for this work was provided by dana-farber cancer institute women’s cancers program, the claudia adams barr foundation and the national cancer institute at the national institutes of health.
to leverage the potential of multi-omics studies, exploratory data analysis methods that provide systematic integration and comparison of multiple layers of omics information are required. we describe multiple co-inertia analysis, an exploratory data analysis method that identifies co-relationships between multiple high dimensional datasets. based on a covariance optimization criterion, mcia simultaneously projects several datasets into the same dimensional space, transforming diverse sets of features onto the same scale, to extract the most variant from each dataset and facilitate biological interpretation and pathway analysis. we demonstrate integration of multiple layers of information using mcia, applied to two typical “omics” research scenarios. the integration of transcriptome and proteome profiles of cells in the nci- cancer cell line panel revealed distinct, complementary features, which together increased the coverage and power of pathway analysis. our analysis highlighted the importance of the leukemia extravasation signaling pathway in leukemia that was not highly ranked in the analysis of any individual dataset. we believe mcia is an attractive method for data integration and visualization of several datasets of multi-omics features observed on the same set of individuals.
 to leverage the potential of multi-omics studies, exploratory data analysis methods that provide systematic integration and comparison of multiple layers of omics information are required. we describe multiple co-inertia analysis, an exploratory data analysis method that identifies co-relationships between multiple high dimensional datasets. based on a covariance optimization criterion, mcia simultaneously projects several datasets into the same dimensional space, transforming diverse sets of features onto the same scale, to extract the most variant from each dataset and facilitate biological interpretation and pathway analysis.  we demonstrate integration of multiple layers of information using mcia, applied to two typical “omics” research scenarios. the integration of transcriptome and proteome profiles of cells in the nci- cancer cell line panel revealed distinct, complementary features, which together increased the coverage and power of pathway analysis. our analysis highlighted the importance of the leukemia extravasation signaling pathway in leukemia that was not highly ranked in the analysis of any individual dataset. secondly, we compared transcriptome profiles of high grade serous ovarian tumors that were obtained, on two different microarray platforms and next generation rna-sequencing, to identify the most informative platform and extract robust biomarkers of molecular subtypes. we discovered that the variance of rna-sequencing data processed using rpkm had greater variance than that with mapsplice and rsem. we provided novel markers highly associated to tumor molecular subtype combined from four data platforms. mcia is implemented and available in the r/bioconductor “omicade4” package.  we believe mcia is an attractive method for data integration and visualization of several datasets of multi-omics features observed on the same set of individuals. the method is not dependent on feature annotation, and thus it can extract important features even when there are not present across all datasets. mcia provides simple graphical representations for the identification of relationships between large datasets. multivariate analysismultiple co-inertiadata integrationomic datavisualization 