BACKGROUND
the constitutive photomorphogenesis   <dig>  or cop <dig> signalosome, is a critical multi-functional protein complex in cells  <cit> .
cop <dig>  originally discovered in arabidopsis thalania as a negative regulator of photomorphogenesis, was later found to be highly conserved across many species  <cit> .
cop <dig> is a 450– <dig> kda holocomplex composed of eight subunits whose official gene symbols are listed in table  <dig>  <cit> .
cop <dig> has catalytic activity to remove nedd <dig>  a small-ubiquitin like modifier, to regulate cullin-ring ubiquitin e <dig> ligases, and thus protein degradation pathways mediated by cullin complexes  <cit> .
cop <dig> regulates phosphorylation of proteins such as tumor suppressor p <dig>  and neddylation essential to ubiquitin-mediated proteolysis of key proteins including tumor suppressor hif-1α .
thus, cop <dig> has been an important focus in cancer genetics .table  <dig> nomenclature of mammalian cop <dig> genes
despite intense studies in the past, the exact essential role of cop <dig> still remains unknown.
this is largely due to cop9’s promiscuous property to interact with a variety of cullin complexes and others, and to affect stabilities of a number of cellular factors.
additionally, subforms of cop <dig> complexes were reported  <cit> .
dubiel et al.
proposed that such flexible conformation enables the broad range of substrates.
among the eight subunits, cops <dig> is the most extensively studied subunit that has numerous critical functions including deneddylation reactions  <cit> .
it is commonly over-expressed in cancers leading to increased deneddylation within cells  <cit> .
other subunits have been shown to possess unique functions.
for example, csn <dig>  is involved in stabilization of p <dig>  <cit> .
over-expression of cops <dig> is linked to chromosome instability  <cit> .
both cops <dig> and cops <dig> knockout are embryonic lethal at e <dig>  and e <dig>  respectively  <cit> .
cops <dig> is frequently over-expressed in cancer  <cit> .
cops <dig> knockout is also embryonic lethal  with impaired growth and differentiation  <cit> .
although these reports convincingly indicate that all the cop <dig> subunits have essential functions in mammalian cells, how these reported functions contribute to cop <dig> activity remains largely unknown.
identification of functional groups with which cop <dig> signalosome may interact in normal physiology may provide an important clue to understand how dysregulation of cop <dig> may facilitate tumorigenesis.
deficiency in even one subunit leads to cop <dig> destabilization.
down-regulation of cops <dig> in mouse embryonic fibroblasts drastically decreased the stability of cops <dig>  cops <dig>  and cops <dig> .
therefore, it is plausible that coordinated expression of cop <dig> genes is needed to maintain necessary amounts of cop <dig> and imbalance may be associated with disease  <cit> .
to elucidate this possibility, a statistically adequate number of normal tissues must be analyzed with precision for their gene expressions.
using cultured cell lines is not appropriate for this investigation, because the cell lines established in vitro are all transformed and thus expressions of the cop <dig> genes are likely altered from cells in normal in vivo conditions.
the cancer genome atlas  has procured specimens from more than  <dig> human tumor tissues from various cancers for comprehensive analyses of tumor tissues and advanced cancer genomics.
importantly, matched normal tissues have been analyzed for about 10% of each cancer type.
using the rnaseq expression data in tcga, we elucidated the possibility that expression of the cop <dig> genes are synchronized in the normal cells.
indeed the cop <dig> gene expressions were in synergy in the normal tissues, but this fine coordination was diminished in the tumor tissues.
furthermore, the analysis led to the unexpected observation that the expressions of the cop <dig> genes were highly synchronized with those involved in mitochondrial biogenesis, particularly those of respiratory chain complex.
methods
nomenclature of cop <dig> genes
although the gene names from csn <dig> to csn <dig> are intuitive and have been accepted as the names for the cop <dig> subunits  <cit> , the tcga uses the official gene symbols.
for the simplicity, this study used the official gene symbols as the names of cop <dig> subunits .
data obtained from tcga
the mrna expression data of normal and tumor tissues were obtained through tcga’s online data portal site .
the cancer tissues chosen for this study are  <dig> cases of head and neck squamous cell carcinoma with  <dig> normal tissues, as well as lung adenocarcinoma .
normalized gene expression results  were collected.
each set of mrna expression data contained  <dig>  gene entries, which were linked to annotations through david functional gene annotation database , yielding  <dig>  analyzable gene lists.
the rnaseq data were linked to information about the chromosome cytoband and gene ontology .
similarly, the de-identified clinical information including anatomical sites and smoking history were retrieved from tcga.
data processing and analysis were carried out using r  and series of command line programs developed with the swift language using xcode ver.
<dig>
generation of correlation coefficient table specific for individual genes
gene expression data of normal and tumor tissues from donors were assembled into a single table file consisting of columns of tissues and rows of genes.
genes that show zero mrna values in more than 5% of the tissues analyzed were excluded unless noted otherwise.
pearson’s correlation coefficient  values between expressions of a target gene  and each of all the other genes were determined, and the genes were then sorted based on the r values.
significance of r values were assessed with p values calculated with cor.test function in r. pair and boxplots were drawn for particular gene sets using r.
expressional correlation of genes closely localized on the chromosomes
generation of cytoband plots
genes ranked within  <dig> highest |r| values were pooled, and the cytoband information was processed into float values by the simple formula as following.
for a gene with cytoband of apb.n or aqb.n, a float value was generated by a - / <dig> or a + / <dig>  respectively.
for example, 8p <dig>  was converted to 8– <dig> /100 =  <dig> , and 8q <dig>  was to  <dig> .
values of  =  was plotted using r. x and y chromosomes were shown in one column as all analyzed genes mapped in x chromosome were also found in y chromosome.
consideration of exact transcription start sites
transcriptional start sites of the cop <dig> genes were retrieved from the ensemble database.
for a particular genea in the rnaseq data, a group of genes  localized within  <dig> mbp of the genea were pooled .
then, an enrichment score  for the gi was calculated based on the correlation coefficient table for genea  <cit> .
prior to this calculation, genea itself  was removed from the table.
separately, genes localized on the same chromosome but outside of the  <dig> mbp of genea were pooled  and used to calculate control mean es and standard deviation  for the nh number of genes that were randomly selected from the go gene pool.
this calculation was repeated for  <dig> times to obtain mean and sd values which were used to generate a p value based on pnorm function with r. to compare these es values with other genes, the above calculation was performed with all available genes in the rnaseq data, and the p-values were used to calculate false discovery rates  using p.adjust function in r.
interpretation of functional gene annotation
the top  <dig> genes with the highest |r| values were pooled with each gene accompanying gene ontology  and kegg pathway information.
over-represented pathways are sorted by the p values based on binomial distribution:  <dig> - pbinom , where n is the actual count occurring in the  <dig> genes, and μ is the average of finding the particular pathways, and calculated with the whole gene set of  <dig> .
derivation of enrichment score  for kegg functional pathways
es scores were calculated by the method described previously  <cit> .
for example, an es for cops <dig> to oxidative phosphorylation pathway was calculated using the cops <dig> correlation coefficient table with gene annotations based on the association with “hsa00190”  to generate phit - pmiss.
es was provided as the maximum value of phit - pmiss as described  <cit> .
a histogram of es values was calculated with randomly chosen nh genes for iteration of  <dig>  and the distribution was used to assess the statistical significance.
ess and corresponding p values were calculated for all genes and then false discovery rates were determined as described above.
RESULTS
a single gene expression analysis
cops <dig>  is responsible for the deneddylation activity of the cop <dig> signalosome, and also known as jab <dig> to activate c-jun .
thus, most cop <dig> studies have focused on the activity of cops <dig>  <cit>  and its influence on biological activities including tumorigenesis.
because of the pivotal role of cops <dig> for the cop <dig> signalosome in mammalian cells, this study primarily focused on cops <dig> to investigate the possibility of coordinated expression of the cop <dig> subunits in the normal and tumor tissues.
our primary interest was the head and neck squamous cell carcinoma , where effects of smoking and influences of major genetic factors, including p <dig>  p <dig>  akt <dig> and k-ras could be addressed due to its established cancer etiology involving these genes.
in addition, hnscc data in tcga included a few distinct anatomical sites including oral cavity, tongue and larynx, and thus effects of tissue sites for the gene expression could be elucidated with the clarified clinical record.
for the simplicity, the official gene symbols will be used in this study instead of more commonly used protein names .
expression of the cops <dig> gene was examined using the data of  <dig> normal head and neck tissues available in tcga.
cops <dig> expression in normal tissues showed mean value of  <dig>  with standard deviation   <dig>  and coefficient of variation   <dig> .
as a comparison, mean coefficient of variation of the  <dig>  gene set was calculated, and found to be  <dig>  .
therefore, the cops <dig> gene expression appeared to be relatively stable among normal oral tissues.
coordinated expression of the cop <dig> genes in normal tissues
expression patterns of cops <dig> and other cop <dig> genes were elucidated by comparing correlation coefficients , using the rna expression data .
high positive correlations were found among the cop <dig> gene expressions.
specifically, relative to cops <dig>  all the other cop <dig> genes except cops7b showed statistically significant positive r values .
gps <dig>  showed a smaller r value than the others, namely cops <dig>  cops <dig>  cops <dig>  cops <dig>  cops7a, and cops <dig>  the lower correlation of gene expression of gps <dig> vs. cops <dig> was found to be a general characteristic through our analyses in this study.
however, the positive correlation between gps <dig> and cops <dig> was still statistically significant .
in addition to the low r value relative to the other cop <dig> genes, the level of cops7b expression was significantly lower compared to cops7a .table  <dig> correlation of expression of cop <dig> genes in normal oral tissues
the rna expression data obtained from tcga were analyzed for pearson’s correlation coefficients.
the number of specimens =  <dig>  p values for cops <dig> vs. other cop <dig> genes were calculated as described based on paired samples algorithm in r
fig.
<dig> pair-wise plot of expressions of cop <dig> genes.
the rna expression results from normal  <dig> oral tissues are plotted.
1: gps <dig>  2: cops <dig>  3: cops <dig>  4: cops <dig>  5: cops <dig>  6: cops <dig>  7a: cops7a, 7b: cops7b, 8: cops8
these results indicated that the expressions of cop <dig> genes were synchronized in human normal oral tissues.
leppert et al.
recently observed that mrna levels of cop <dig> genes were regulated by let- <dig> mirna  <cit> .
to probe the importance of the regulatory mechanism by mirna, the tcga mirna dataset was used to determine the gene expression correlation of cop <dig> genes and the let- <dig> mirna species .
there were statistically significant negative correlations between let- <dig> isoforms  and cop <dig> genes , which was consistent with the previous work.
loss of synergistic expressions of the cop <dig> genes in cancer tissues
to determine whether the synergistic expressions of the cop <dig> genes are maintained in the cancer tissues,  <dig> tumor tissues matched to the normal cases were pooled , and determined their correlation coefficients .
high positive r values were observed among cop <dig> genes in normal tissues as expected.
in contrast, the correlations of expressions among cop <dig> genes diminished in the matched tumor tissues .
moreover, the negative effect of let- <dig> mirna observed in the normal tissue dataset was diminished in the matched tumor tissues .
these results demonstrated that the coordinated expression of cop <dig> genes was diminished in these tumor tissues.table  <dig> loss of coordinated cop <dig> gene expression in hnscc
normal and matched hnscc tissues  were compared for their correlation coefficient values among the cop <dig> genes
validation of the results with subgroups
hnscc is an age-related disease.
the majority of specimens in tcga are from elderly patients , and less than 4% of the patients are younger than 40 years old.
thus, the loss of the coordination of the cop <dig> gene expressions may be due to the age-related degeneration of the coordination.
however, we analyzed the cop <dig> expression in the  <dig> tumor tissues that were 42 years old or younger to the  <dig> tumor tissues that were 80 years old or older .
with this analysis, we found similar loss of the coordinated expressions in specimens from both the young and elderly patients.
in addition, there were no significant differences of the expression levels of each cop <dig> gene between the young and aged groups .
these results support that the age is not a factor causing the low coordination of the cop <dig> genes in the tumor specimens.
the oral tissues in tcga were taken from different anatomical sites .
the synergistic expression of cop <dig> subunit genes may vary depending on the exact site of oral tissues.
to assess this possibility, the data set was arranged into subgroups based on the sites.
the major groups were the oral cavity  and the oral tongue .
the data sets of the two subgroups were used to determine correlation coefficients of the cop <dig> genes.
these values were similar to those of the entire dataset , indicating that the synergistic expression of the cop <dig> genes are maintained regardless of the sites of the oral tissues.
similarly, there appeared to be no clear effects of smoking, because analyzing the normal tissues from eleven current smokers identified a highly similar synergistic rna expressions of cop <dig> genes .
chromosome mapping of genes of high synergy with the cops <dig> expression
to elucidate the tumor-associated disintegration of the cop <dig> gene expression in detail, we focused on the cops <dig> gene.
we generated a “correlation table” of the cops <dig>  in which expressional correlation coefficients relative to the cops <dig> gene  were calculated for all the analyzable genes .
by cross-referring the results to gene ontology provided with david, we noticed that genes with high r
cops <dig> of the tumor tissues were mapped in the same chromosomal location as the cops <dig> gene .
we thus speculated that the loss of the coordination of the cop <dig> expressions in the tumor tissues was not due to unpredictably disorganized expressions of the cop <dig> subunit genes, but rather there was a dominant effect of the physical locations of the genes on their expressions.
to illustrate this possibility, chromosomal locations of  <dig> genes with top | r
cops5| values were plotted on chromosomes based on the cytoband information .
the genes with high |r
cops5| in the normal tissues were distributed through all the chromosomes; among the top  <dig> genes, only seven genes were mapped at 8q, near the cops <dig> chromosomal location .
in contrast, genes with high |r
cops5| in the tumor tissues were concentrated near cops5’s chromosomal region ;  <dig> genes out of the top  <dig> genes were bound to 8q.fig.
<dig> chromosomal mapping of genes with highly synchronized expression with the cops <dig> gene.
the  <dig> highest ranked genes based on |rcops5| were pooled, and mapped on to chromosomal locations based on their cytoband information.
in the top panel, a cartoon denotes the chromosome  <dig>  where the cops <dig> gene is located .
the major grid marks on the x-axis  depict location of centromeres dividing p and q arms at the left and the right sides, respectively.
x and y chromosomes are combined, because all genes analyzed in the figure are mapped in both x and y
we sought to assess the chromosomal cis-effect on the cops <dig> gene expression more rigorously.
to this end, based on the method by subramanian et al.
we determined the enrichment score  based on the r
cops <dig> for the groups of genes that are localized within  <dig> mega base pairs  from the cops <dig>  and obtained a p-value based on the null hypothesis that these es values were not different from that calculated for genes outside  <dig> mbp from the cops <dig>  <cit> .
we have also determined such es values for all the genes in the chromosome  <dig> to determine the false discovery rate .
the es value for the cops <dig> with the tumor tissues  was significantly higher than that with the normal tissues .
the procedure was applied for the other cop <dig> genes, and revealed significant links to the chromosomal cis-effect in the case of the tumor tissues .
therefore, we conclude that the cop <dig> expression in the tumors is disintegrated because the chromosomal cis-effect becomes a dominant force to express each subunit gene.
furthermore, these results indicated that analysis of the cop <dig> expressions in the normal tissues rather than tumor tissues is pivotal for understanding the cop9’s role in vivo, and such an analysis is possible using the rnaseq data in the normal tissues available from tcga.
pathway analysis based on kegg and go revealed association of cops <dig> expression with mitochondrial pathways
genes that belong to particular cell pathways may be over-represented in the cops <dig> correlation coefficient table.
identifying such pathways should provide important clues as to cop9’s roles under the normal physiology.
based on the results above, such functional association should be elucidated with the data derived from normal tissues, and the large pool of high quality gene expression data in the tcga database prompted us to investigate this possibility.
genes with the top  <dig> |r
cops5| values in the normal tissues were pooled, and functional pathways available through go and kegg database were associated to each gene .
the go and kegg pathways were then ranked based on the over-representation, and the  <dig> most over-represented pathways are shown in table  <dig> .
unexpectedly, the pathways related to mitochondrial functions and components were highly over-represented.
for example,  <dig> genes belonging to the oxidative phosphorylation kegg pathway  were among the genes with high expression correlation with cops <dig> out of  <dig> total genes .
a number of go pathways involving mitochondrial functions were also found among the over-represented pathways.
these included respiratory chain , mitochondrial electron transport , nadh dehydrogenase activity , mitochondrial respiratory chain complex i .
in fact, all of the  <dig> most over-represented pathways are directly related to major mitochondrial and energy derivation functions.
at the time this study was carried out, gene ontology database had not yet linked cops <dig> to any of mitochondrial-related pathways.
expression values were examined in pair-plots for cops <dig> and mitochondrial genes of the five highest r
cops <dig>  ndufb <dig>  mrps <dig>  atp5f <dig>  tmem126a, ndufb <dig> .
the results illustrate high integrity of the data and the coordination of the expression of the cops <dig> with the mitochondria related genes.table  <dig> association of cops <dig> with mitochondrial pathways in normal oral tissues
found: number of genes in a particular go or kegg pathway found in the top  <dig> r
cops <dig> list; expected: number of expected genes to yield p <  <dig>  appearing in  <dig> genes based on the average appearance per gene ; μ values were calculated by the total appearance of individual pathways in the entire gene list of  <dig> .
p values were calculated using binomial distribution in r. expected occurrences in a list of  <dig> genes to provide p <  <dig>  are listed
fig.
<dig> pair-wise plot for rna expression  <dig> mitochondrial genes.
the five genes with the highest r
cops <dig> in the normal oral tissues are plotted.
1: cops <dig>  2: ndufb <dig>  3: mrps <dig>  4: atp5f <dig>  5: tmem126a, 6: ndufb3
dependency on the anatomical sites and smoking history were tested with the same subgroups as described above .
the mrna expression data from the two major anatomical sites, oral cavity and oral tongue, resulted in almost identical over-representation of the mitochondrial pathways with cops <dig>  the same conclusion was made with the subgroup containing only the current smokers.
alteration of cops <dig> expression coordination in cancer tissues
the gene expression coordination between the cops <dig> gene and the mitochondrial pathways was elucidated with tumor tissues .
although a few go pathways related to mitochondrial functions were over-represented, major pathways associated with the cops <dig> expression were those involved in ribosome processing and rna binding as well as membrane lumen pathways.table  <dig> loss of coordinated expression of cops <dig> with mitochondria related genes in the tumor tissues
the analysis was carried out in the same way as for table 5
the possibility that tissues from a particular anatomical site with high dysregulation caused the bias was addressed.
oral cavity, tongue, and larynx are the major three sites that together contribute to about 80% of the donated head and neck cancer tissues.
cops5-associated go and kegg pathways were determined for the tissues from the particular sites, and confirmed that correlation of cops <dig> expression with mitochondria associated genes were not as predominant as in the normal tissues.
to further support the finding, enrichment score  was calculated for the genes in the oxidative phosphorylation kegg pathway  using the rnaseq data of  <dig> normal head and neck tissues  <cit>  .
the obtained es value,  <dig> , was much higher than the average es value based on randomly assigned gene pools .fig.
<dig> enrichment score for kegg oxidative phosphorylation pathway.
a the enrichment score was calculated for the oxidative phosphorylation kegg pathway .
genes are aligned on x-axis based on their correlation coefficients to the cops <dig> expression.
the hsa <dig> pathway contained  <dig> genes in the entire list.
b  <dig> randomly selected genes were used to calculate es using the r
cops <dig> gene list as in a. the histogram was generated by iteration of  <dig> times.
the arrow indicates the es for hsa <dig>
in addition, es values corresponding to all the other kegg pathways were calculated for the  <dig> matched normal and tumor tissues .
the es value of the oxidative phosphorylation pathway  found to be among the highest.
the ubiquinone biosynthesis pathway generated the highest es value, which is an essential factor in the respiratory chain and requires phenylalanine  for its synthesis.
thus, the enrichment of genes involved in mitochondria-related pathway based on the cops <dig> correlation was exceptionally high.table  <dig> enrichment scores for the cops <dig> gene on kegg pathways
all kegg pathways found in the entire gene list of the tcga rna expression database  were individually scanned through the r
cops <dig> alignment, and the maximum  as well as minimum  were calculated as previously described  <cit> .
the kegg pathways that generated the  <dig> highest es values are shown.
nh: number of genes belonging to a kegg pathway in the entire list.
a p value for a particular kegg pathway with an nh  was determined using mean and standard deviation obtained by the following method.
nk genes were randomly selected from the r
cops <dig> list and calculated an es value, which was reiterated for  <dig> times to generate the mean and standard deviation of es
to elucidate the significance of the cop <dig> coordination with oxidative phosphorylation compared to other genes, we calculated es values relative to the hsa <dig> kegg pathway for the entire gene list in the rnaseq data of the matched normal and tumor tissues to calculate fdr .
except for the cops7b gene, all the cop <dig> subunit genes showed highly significant correlations with the genes in the oxidative phosphorylation pathway.
not surprisingly, most of the genes in the hsa <dig> oxidative phosphorylation pathway resulted in the high es/low fdr values, validating our approach to detect a linkage of a gene to the oxidative phosphorylation.
remarkably, the cop <dig> genes  were ranked comparably with the hsa <dig> genes.
the p as well as fdr values were increased for all the cop <dig> genes  with the tumors, indicating significant decrease in the association of the cop <dig> expression to the oxidative phosphorylation related genes.
taken all together, we conclude that there is a high synergy between the expression of the cop <dig> and mitochondria functions, which implies that the normal physiological function of the cop <dig> may be inseparable from the mitochondrial activity.
analysis of gene expression of other cancers in the tcga database
finally, the question was addressed whether the association of the expressions of cop <dig> genes with those of mitochondrial genes was unique in head and neck tissues.
to probe this, the rna expression data of tcga for the lung squamous cell carcinoma tissues were analyzed.
tumors originating from the lung upper airway epithelia share similarity with hnscc in its etiology and its involvement of smoking as an important risk factor.
first, the synergistic expression pattern of cop <dig> genes in the normal tissues  was examined, and found to be highly synergistic with one another .
exceptions were gps <dig> and cops7b, which showed lower correlation coefficients to other cop <dig> factors.
the high positive correlations were again lost in the tumor tissues.
these results showed a remarkable coincidence with those of head and neck normal and tumor tissues.table  <dig> synchronized expression of the cop <dig> genes in the matched tissues of normal and lung squamous cell carcinoma
<dig>
<dig>
<dig>
<dig>
<dig>
<dig>
the analysis was carried out in the same way as for table  <dig>  values in bold letters indicate significantly high cosp <dig> correlation with other cop <dig> components
the functional pathway analysis with the entire gene list resulted in the high correlation of the cops <dig> expression with the mitochondrial genes .
unlike the case of head and neck squamous cell carcinoma, however, the coordinated expression of cop <dig> and the mitochondrial genes was maintained well in the lung squamous cell carcinoma tissues .
discussion
cancer genetics and analyses of the tcga database mainly focus on the loss of genomic integrity.
mutations and single nucleotide polymorphism in dna show the high rate of permanent alteration of the genome, epigenetic modification of the genome, rna and mirna expression alterations.
tcga data sets are available publicly for these types of analyses.
the database is highly reliable and allowed us to examine cancer genomics with more than  <dig> tumor specimens  <cit> .
in addition to the tumor specimens, data from about  <dig> to  <dig> normal tissues are made available for most types of cancers with a notable exception for glioblastoma.
although the data sizes are smaller than those for the cancer tissues, they are arguably one of the most comprehensive databases for genomics of human specimens of non-disease status with pathological and clinical scrutiny, which should provide valuable information for elucidating functions of a gene under the normal physiology.
we hypothesized that a synergistic expression of genes may be observed in the normal tissues, when expressions of a group of genes need to be coordinated for their proper cellular activities.
cop <dig> signalosome is an ideal subject to test this possibility, because all eight subunits are essential for the proper cop <dig> function; deficiency in one of the cop <dig> genes decreased stability of the other subunits and caused lethality  <cit> .
in this study the cop <dig> signalosome genes showed highly synergistic expressions in the normal head and neck tissues, except for gps <dig> and cops7b.
remarkably, these characteristics, including the lower r values for gps <dig> and cops7b, were maintained well even in the lung and breast normal tissues .
the coordinated expression of the cop <dig> genes decreased in the tumor tissues.
therefore, this coordination could not have been detected by analyzing tumor tissues.
likewise, because all cells cultured in vitro are under non-physiological conditions and thus transformed at some extent, it would have been difficult to obtain convincing evidence for the expression correlation of the cop <dig> subunit genes based on in vitro cell biological studies alone.
it is interesting to observe the lowest coordination of gps <dig> expression compared to the other cop <dig> components .
besides being a subunit of cop <dig>  gps <dig> is known for its function in suppression of g-protein and mitogen-activated signal transduction in mammalian cells  <cit> .
this cop9-independent function quite likely requires a regulation of the gps <dig> gene separate from the other cop <dig> genes.
cops7b not only had the lowest r value, but also significantly lower expression compared to cops7a .
while cops7a and cops7b genes are encoded in different human chromosomes, chromosome  <dig> and  <dig> respectively, there is a high homology between cops7a and cops7b in the amino acid sequence .
thus, cops7b may have only a minor functional role for cop <dig>  which may be compensated by cops7a.
to our knowledge, using the exact physical map of the transcription start sites available from the ensemble genome database has not been applied previously to elucidate effects of the relative locations of genes on the gene expression correlations.
in this study, we used this approach to find the chromosomal locations of the genes as a dominant factor causing the loss of the coordinated expressions of the cop <dig> genes.
much fewer number of genes were under the control of this cis-effect in case of normal tissues, indicating that most genes in normal tissues are regulated independently of their chromosomal locations and perhaps more based on their functional groups, as supported by this study.
therefore, important implications derived from this analysis are that 1) analyzing tumor tissue data may not identify a normal function of a protein family such as cop <dig>  and 2) the qualities of the normal tissue data available in tcga are high enough to help identifying unknown expressional linkages with functional pathways.
although the primary focus of this study was cop <dig> expression profiling, we expanded the analysis of the influence of the chromosomal location for all the genes that were analyzable in the rnaseq data, and the data for cancer genes  are shown in additional file 7: table s <dig>  it may be interesting to elucidate whether these proto-oncogenes have tight association with the chromosomal locations in other types of tumor tissues available in tcga.
the observation that the expression of the cops <dig> gene  was synchronized with those involved in mitochondrial function was unexpected.
the correlation coefficients of the cops <dig> and top ranked genes such as ndufb <dig> and mrps <dig>  were remarkably high with statistical significance.
moreover, not only a handful of genes involved in mitochondria but also the total pathway analysis revealed that mitochondria-related cellular functions and pathways were predominantly over-represented in the cops <dig> correlation coefficient table.
although there have not been reports that show direct function of cop <dig> for mitochondria, correlation analyses with the tp <dig> of which functions have been extensively studied linked it with known cellular functions including cell cycle and dna damage response pathways .
the result validates the approach in this study, and support the prediction that the cop <dig> signalosome has a novel functional link to mitochondrial activities.
few studies have linked cop <dig> to mitochondria-related cellular activities.
jab1/cops <dig> was previously identified as an apoptosis inducing factor through modulation of mitochondrial bh3-domain containing protein bclgs  <cit> .
a study found that in neurospora crassa, three subunits of the cop <dig> were required to recover from deficiency of the alternative oxidase function, which is important for mitochondrial functions  <cit> .
however, no follow up studies have been carried out for mammalian cells.
on the other hand, pathways involved in ribosomal functions were associated with the cops <dig> expression in the tumor tissues.
the significance of this link is currently unclear.
however, the cop <dig> structure is homologous to translation initiation factor  <dig> , another high-ordered complex in the cells  <cit> , and our speculation is that cops <dig> may have an overlapping function in the protein synthesis in transformed cells.
the quantity and quality of the genomics database will only be improved.
the present study underscores the importance and usefulness of examining the genomics of normal tissues.
the results would provide novel functional links of a gene, which was not possible by analyzing the data set of pathophysiological specimens.
this approach can reasonably be applied to other genes to identify pathways that have not been known due to scarce information for genomics information of normal tissues.
CONCLUSIONS
the present study investigated cops <dig> gene expression using the tcga database.
the analysis revealed the highly synergistic mrna expression pattern among cop <dig> genes in normal oral tissues.
the coordination was abrogated in the tumor tissues, and the regulatory mechanism was taken over by a cis-acting gene regulation.
further analysis revealed an unexpected expression correlation between the cops <dig> gene and a number of mitochondria-related genes, postulating the possible functional role of the cop <dig> signalosome for mitochondrial activities.
the quality of the rna expression data available in tcga is high enough to allow us gene expression analysis in normal tissues to identify gene functions under normal physiology, which may not be recognized by the data obtained from tumor tissues.
additional files
for analysis of rna expression of normal and cancer tissues reveal high correlation of cop <dig> gene expression with respiratory chain complex components by christina a. wicker and tadahide izumi.additional file 1: figure s <dig>  expressions of cop <dig> genes in normal oral tissues.
the rna expression levels  of cop <dig> genes from  <dig> normal oral tissues are shown.
horizontal lines in the boxes: median values.
top and bottom whiskers denote the first and the third quartile .
outliers  and < q <dig> -  <dig>  * ) are shown as individual dots.
figure s <dig>  pair-wise plot for expressions of cop <dig> genes in the matched lung tissues of normal and squamous cell carcinoma.
the plots were generated as described in fig.
<dig>  except for the data set of the normal  and their matched lung squamous cell carcinoma  from tcga was used.
the number of samples =  <dig>
additional file 2: table s <dig>  the let- <dig> expression profilesof the matched head and neck normal and tumor tissues.
correlation coefficients for individual let- <dig> isoforms relative to cop <dig> genes are shown.
asterisks: p <  <dig> .
table s <dig>  to support the result of table  <dig>   <dig> additional normal head and neck tissues were included to calculate correlation of expressions among the cop <dig> genes.
correlation coefficients of cop <dig> genes versus cops <dig>  were determined from randomly selected  <dig> tumor specimens from the entire pool of  <dig>  the calculations were repeated for  <dig>  times to obtain their mean and standard deviation  values, which were used to compare the difference from rcops <dig> of normal tissues.
table s <dig>  rnaseq data from young  and elderly patients  were used to calculate the correlation coefficients among the cop <dig> genes.
table s <dig>   <dig> normal oral cavity tissues were pooled based on the tcga clinical data, and the correlation coefficients and their p values for cops <dig> were calculated as described in the main article.
table s <dig>  data based on the  <dig> normal tongue tissues were pooled based on the clinical data in the tcga database, and the correlation coefficients and their p values for cops <dig> were calculated.
table s <dig>   <dig> normal tissues from smokers were pooled based on the clinical data in the tcga database, and the correlation coefficients and their p values for cops <dig> were calculated.
table s <dig>  nh: number of genes found within  <dig> mbp of each cop <dig> gene; mean, sd, p and q values were determined as described in materials and methods.
table s <dig>  es for cop <dig> genes for all kegg pathways.
table s <dig>  es for cop <dig> genes for the oxidative phosphorylation pathway.
additional file 3: table s <dig>  the cops <dig> correlation coefficient table list includes  <dig> highest |r| values relative to the expression of cops <dig> in the normal and tumor tissues of head and neck.
additional file 4: table s <dig>  go and kegg pathways that are over-represented in the normal and tumor oral tissues based on the cops <dig> correlation coefficient table.
additional file 5: table s <dig>  go and kegg pathways over-represented by cops <dig> expression correlation in the two major anatomical sites in the head and neck .
additional file 6: table s <dig>  go and kegg pathways over-represented in the normal and tumor lung tissues based on the cops <dig> correlation coefficient table.
additional file 7: table s <dig>  genes under the chromosomal cis-effect in the hnscc tissues
additional file 8: table s <dig>  kegg pathway analysis for the tp <dig> gene
funding
this work was supported by the nih  training grant t32es <dig> and nci grants ca <dig>  cancer center supporting grant p30ca <dig>  and ctsa ul1tr <dig>
availability of data and materials
all the data discussed in this study are available as in main figures or in the supplemental data.
authors’ contributions
ti conceived the project idea and wrote necessary programs.
ti and caw performed the data analysis.
ti and caw wrote and edited the manuscript, and approved the final manuscript.
competing interests
the authors declare that they have no competing interests.
consent for publication
not applicable.
ethics approval and consent to participate
not applicable.

 the constitutive photomorphogenesis or cop signalosome, is a critical multi-functional protein complex in cells. cop originally discovered in arabidopsis thalania as a negative regulator of photomorphogenesis, was later found to be highly conserved across many species. cop is a 450– kda holocomplex composed of eight subunits whose official gene symbols are listed in table . cop has catalytic activity to remove nedd a small-ubiquitin like modifier, to regulate cullin-ring ubiquitin e ligases, and thus protein degradation pathways mediated by cullin complexes. cop regulates phosphorylation of proteins such as tumor suppressor p and neddylation essential to ubiquitin-mediated proteolysis of key proteins including tumor suppressor hif-1α. thus, cop has been an important focus in cancer genetics.table nomenclature of mammalian cop genes despite intense studies in the past, the exact essential role of cop still remains unknown. this is largely due to cop9’s promiscuous property to interact with a variety of cullin complexes and others, and to affect stabilities of a number of cellular factors. additionally, subforms of cop complexes were reported. dubiel et al. proposed that such flexible conformation enables the broad range of substrates. among the eight subunits, cops is the most extensively studied subunit that has numerous critical functions including deneddylation reactions. it is commonly over-expressed in cancers leading to increased deneddylation within cells. other subunits have been shown to possess unique functions. for example, csn is involved in stabilization of p. over-expression of cops is linked to chromosome instability. both cops and cops knockout are embryonic lethal at e and e respectively. cops is frequently over-expressed in cancer. cops knockout is also embryonic lethal with impaired growth and differentiation. although these reports convincingly indicate that all the cop subunits have essential functions in mammalian cells, how these reported functions contribute to cop activity remains largely unknown. identification of functional groups with which cop signalosome may interact in normal physiology may provide an important clue to understand how dysregulation of cop may facilitate tumorigenesis. deficiency in even one subunit leads to cop destabilization. down-regulation of cops in mouse embryonic fibroblasts drastically decreased the stability of cops cops and cops. therefore, it is plausible that coordinated expression of cop genes is needed to maintain necessary amounts of cop and imbalance may be associated with disease. to elucidate this possibility, a statistically adequate number of normal tissues must be analyzed with precision for their gene expressions. using cultured cell lines is not appropriate for this investigation, because the cell lines established in vitro are all transformed and thus expressions of the cop genes are likely altered from cells in normal in vivo conditions. the cancer genome atlas has procured specimens from more than human tumor tissues from various cancers for comprehensive analyses of tumor tissues and advanced cancer genomics. importantly, matched normal tissues have been analyzed for about 10% of each cancer type. using the rnaseq expression data in tcga, we elucidated the possibility that expression of the cop genes are synchronized in the normal cells. indeed the cop gene expressions were in synergy in the normal tissues, but this fine coordination was diminished in the tumor tissues. furthermore, the analysis led to the unexpected observation that the expressions of the cop genes were highly synchronized with those involved in mitochondrial biogenesis, particularly those of respiratory chain complex. methods nomenclature of cop genes although the gene names from csn to csn are intuitive and have been accepted as the names for the cop subunits, the tcga uses the official gene symbols. for the simplicity, this study used the official gene symbols as the names of cop subunits. data obtained from tcga the mrna expression data of normal and tumor tissues were obtained through tcga’s online data portal site. the cancer tissues chosen for this study are cases of head and neck squamous cell carcinoma with normal tissues, as well as lung adenocarcinoma. normalized gene expression results were collected. each set of mrna expression data contained gene entries, which were linked to annotations through david functional gene annotation database, yielding analyzable gene lists. the rnaseq data were linked to information about the chromosome cytoband and gene ontology. similarly, the de-identified clinical information including anatomical sites and smoking history were retrieved from tcga. data processing and analysis were carried out using r and series of command line programs developed with the swift language using xcode ver. generation of correlation coefficient table specific for individual genes gene expression data of normal and tumor tissues from donors were assembled into a single table file consisting of columns of tissues and rows of genes. genes that show zero mrna values in more than 5% of the tissues analyzed were excluded unless noted otherwise. pearson’s correlation coefficient values between expressions of a target gene and each of all the other genes were determined, and the genes were then sorted based on the r values. significance of r values were assessed with p values calculated with cor.test function in r. pair and boxplots were drawn for particular gene sets using r. expressional correlation of genes closely localized on the chromosomes generation of cytoband plots genes ranked within highest |r| values were pooled, and the cytoband information was processed into float values by the simple formula as following. for a gene with cytoband of apb.n or aqb.n, a float value was generated by a - / or a + / respectively. for example, 8p was converted to 8– /100 = , and 8q was to. values of  =  was plotted using r. x and y chromosomes were shown in one column as all analyzed genes mapped in x chromosome were also found in y chromosome. consideration of exact transcription start sites transcriptional start sites of the cop genes were retrieved from the ensemble database. for a particular genea in the rnaseq data, a group of genes localized within mbp of the genea were pooled. then, an enrichment score for the gi was calculated based on the correlation coefficient table for genea. prior to this calculation, genea itself was removed from the table. separately, genes localized on the same chromosome but outside of the mbp of genea were pooled and used to calculate control mean es and standard deviation for the nh number of genes that were randomly selected from the go gene pool. this calculation was repeated for times to obtain mean and sd values which were used to generate a p value based on pnorm function with r. to compare these es values with other genes, the above calculation was performed with all available genes in the rnaseq data, and the p-values were used to calculate false discovery rates using p.adjust function in r. interpretation of functional gene annotation the top genes with the highest |r| values were pooled with each gene accompanying gene ontology and kegg pathway information. over-represented pathways are sorted by the p values based on binomial distribution: - pbinom, where n is the actual count occurring in the genes, and μ is the average of finding the particular pathways, and calculated with the whole gene set of. derivation of enrichment score for kegg functional pathways es scores were calculated by the method described previously. for example, an es for cops to oxidative phosphorylation pathway was calculated using the cops correlation coefficient table with gene annotations based on the association with “hsa00190” to generate phit - pmiss. es was provided as the maximum value of phit - pmiss as described. a histogram of es values was calculated with randomly chosen nh genes for iteration of and the distribution was used to assess the statistical significance. ess and corresponding p values were calculated for all genes and then false discovery rates were determined as described above.  a single gene expression analysis cops is responsible for the deneddylation activity of the cop signalosome, and also known as jab to activate c-jun. thus, most cop studies have focused on the activity of cops and its influence on biological activities including tumorigenesis. because of the pivotal role of cops for the cop signalosome in mammalian cells, this study primarily focused on cops to investigate the possibility of coordinated expression of the cop subunits in the normal and tumor tissues. our primary interest was the head and neck squamous cell carcinoma, where effects of smoking and influences of major genetic factors, including p p akt and k-ras could be addressed due to its established cancer etiology involving these genes. in addition, hnscc data in tcga included a few distinct anatomical sites including oral cavity, tongue and larynx, and thus effects of tissue sites for the gene expression could be elucidated with the clarified clinical record. for the simplicity, the official gene symbols will be used in this study instead of more commonly used protein names. expression of the cops gene was examined using the data of normal head and neck tissues available in tcga. cops expression in normal tissues showed mean value of with standard deviation and coefficient of variation. as a comparison, mean coefficient of variation of the gene set was calculated, and found to be. therefore, the cops gene expression appeared to be relatively stable among normal oral tissues. coordinated expression of the cop genes in normal tissues expression patterns of cops and other cop genes were elucidated by comparing correlation coefficients, using the rna expression data. high positive correlations were found among the cop gene expressions. specifically, relative to cops all the other cop genes except cops7b showed statistically significant positive r values. gps showed a smaller r value than the others, namely cops cops cops cops cops7a, and cops the lower correlation of gene expression of gps vs. cops was found to be a general characteristic through our analyses in this study. however, the positive correlation between gps and cops was still statistically significant. in addition to the low r value relative to the other cop genes, the level of cops7b expression was significantly lower compared to cops7a.table correlation of expression of cop genes in normal oral tissues the rna expression data obtained from tcga were analyzed for pearson’s correlation coefficients. the number of specimens =  p values for cops vs. other cop genes were calculated as described based on paired samples algorithm in r fig. pair-wise plot of expressions of cop genes. the rna expression results from normal oral tissues are plotted. 1: gps 2: cops 3: cops 4: cops 5: cops 6: cops 7a: cops7a, 7b: cops7b, 8: cops8 these results indicated that the expressions of cop genes were synchronized in human normal oral tissues. leppert et al. recently observed that mrna levels of cop genes were regulated by let- mirna. to probe the importance of the regulatory mechanism by mirna, the tcga mirna dataset was used to determine the gene expression correlation of cop genes and the let- mirna species. there were statistically significant negative correlations between let- isoforms and cop genes, which was consistent with the previous work. loss of synergistic expressions of the cop genes in cancer tissues to determine whether the synergistic expressions of the cop genes are maintained in the cancer tissues, tumor tissues matched to the normal cases were pooled, and determined their correlation coefficients. high positive r values were observed among cop genes in normal tissues as expected. in contrast, the correlations of expressions among cop genes diminished in the matched tumor tissues. moreover, the negative effect of let- mirna observed in the normal tissue dataset was diminished in the matched tumor tissues. these results demonstrated that the coordinated expression of cop genes was diminished in these tumor tissues.table loss of coordinated cop gene expression in hnscc normal and matched hnscc tissues were compared for their correlation coefficient values among the cop genes validation of the results with subgroups hnscc is an age-related disease. the majority of specimens in tcga are from elderly patients, and less than 4% of the patients are younger than 40 years old. thus, the loss of the coordination of the cop gene expressions may be due to the age-related degeneration of the coordination. however, we analyzed the cop expression in the tumor tissues that were 42 years old or younger to the tumor tissues that were 80 years old or older. with this analysis, we found similar loss of the coordinated expressions in specimens from both the young and elderly patients. in addition, there were no significant differences of the expression levels of each cop gene between the young and aged groups. these results support that the age is not a factor causing the low coordination of the cop genes in the tumor specimens. the oral tissues in tcga were taken from different anatomical sites. the synergistic expression of cop subunit genes may vary depending on the exact site of oral tissues. to assess this possibility, the data set was arranged into subgroups based on the sites. the major groups were the oral cavity and the oral tongue. the data sets of the two subgroups were used to determine correlation coefficients of the cop genes. these values were similar to those of the entire dataset, indicating that the synergistic expression of the cop genes are maintained regardless of the sites of the oral tissues. similarly, there appeared to be no clear effects of smoking, because analyzing the normal tissues from eleven current smokers identified a highly similar synergistic rna expressions of cop genes. chromosome mapping of genes of high synergy with the cops expression to elucidate the tumor-associated disintegration of the cop gene expression in detail, we focused on the cops gene. we generated a “correlation table” of the cops in which expressional correlation coefficients relative to the cops gene were calculated for all the analyzable genes. by cross-referring the results to gene ontology provided with david, we noticed that genes with high r cops of the tumor tissues were mapped in the same chromosomal location as the cops gene. we thus speculated that the loss of the coordination of the cop expressions in the tumor tissues was not due to unpredictably disorganized expressions of the cop subunit genes, but rather there was a dominant effect of the physical locations of the genes on their expressions. to illustrate this possibility, chromosomal locations of genes with top | r cops5| values were plotted on chromosomes based on the cytoband information. the genes with high |r cops5| in the normal tissues were distributed through all the chromosomes; among the top genes, only seven genes were mapped at 8q, near the cops chromosomal location. in contrast, genes with high |r cops5| in the tumor tissues were concentrated near cops5’s chromosomal region; genes out of the top genes were bound to 8q.fig. chromosomal mapping of genes with highly synchronized expression with the cops gene. the highest ranked genes based on |rcops5| were pooled, and mapped on to chromosomal locations based on their cytoband information. in the top panel, a cartoon denotes the chromosome where the cops gene is located. the major grid marks on the x-axis depict location of centromeres dividing p and q arms at the left and the right sides, respectively. x and y chromosomes are combined, because all genes analyzed in the figure are mapped in both x and y we sought to assess the chromosomal cis-effect on the cops gene expression more rigorously. to this end, based on the method by subramanian et al. we determined the enrichment score based on the r cops for the groups of genes that are localized within mega base pairs from the cops and obtained a p-value based on the null hypothesis that these es values were not different from that calculated for genes outside mbp from the cops. we have also determined such es values for all the genes in the chromosome to determine the false discovery rate. the es value for the cops with the tumor tissues was significantly higher than that with the normal tissues. the procedure was applied for the other cop genes, and revealed significant links to the chromosomal cis-effect in the case of the tumor tissues. therefore, we conclude that the cop expression in the tumors is disintegrated because the chromosomal cis-effect becomes a dominant force to express each subunit gene. furthermore, these results indicated that analysis of the cop expressions in the normal tissues rather than tumor tissues is pivotal for understanding the cop9’s role in vivo, and such an analysis is possible using the rnaseq data in the normal tissues available from tcga. pathway analysis based on kegg and go revealed association of cops expression with mitochondrial pathways genes that belong to particular cell pathways may be over-represented in the cops correlation coefficient table. identifying such pathways should provide important clues as to cop9’s roles under the normal physiology. based on the results above, such functional association should be elucidated with the data derived from normal tissues, and the large pool of high quality gene expression data in the tcga database prompted us to investigate this possibility. genes with the top |r cops5| values in the normal tissues were pooled, and functional pathways available through go and kegg database were associated to each gene. the go and kegg pathways were then ranked based on the over-representation, and the most over-represented pathways are shown in table . unexpectedly, the pathways related to mitochondrial functions and components were highly over-represented. for example, genes belonging to the oxidative phosphorylation kegg pathway were among the genes with high expression correlation with cops out of total genes. a number of go pathways involving mitochondrial functions were also found among the over-represented pathways. these included respiratory chain, mitochondrial electron transport, nadh dehydrogenase activity, mitochondrial respiratory chain complex i. in fact, all of the most over-represented pathways are directly related to major mitochondrial and energy derivation functions. at the time this study was carried out, gene ontology database had not yet linked cops to any of mitochondrial-related pathways. expression values were examined in pair-plots for cops and mitochondrial genes of the five highest r cops ndufb mrps atp5f tmem126a, ndufb. the results illustrate high integrity of the data and the coordination of the expression of the cops with the mitochondria related genes.table association of cops with mitochondrial pathways in normal oral tissues found: number of genes in a particular go or kegg pathway found in the top r cops list; expected: number of expected genes to yield p <  appearing in genes based on the average appearance per gene; μ values were calculated by the total appearance of individual pathways in the entire gene list of. p values were calculated using binomial distribution in r. expected occurrences in a list of genes to provide p <  are listed fig. pair-wise plot for rna expression mitochondrial genes. the five genes with the highest r cops in the normal oral tissues are plotted. 1: cops 2: ndufb 3: mrps 4: atp5f 5: tmem126a, 6: ndufb3 dependency on the anatomical sites and smoking history were tested with the same subgroups as described above. the mrna expression data from the two major anatomical sites, oral cavity and oral tongue, resulted in almost identical over-representation of the mitochondrial pathways with cops the same conclusion was made with the subgroup containing only the current smokers. alteration of cops expression coordination in cancer tissues the gene expression coordination between the cops gene and the mitochondrial pathways was elucidated with tumor tissues. although a few go pathways related to mitochondrial functions were over-represented, major pathways associated with the cops expression were those involved in ribosome processing and rna binding as well as membrane lumen pathways.table loss of coordinated expression of cops with mitochondria related genes in the tumor tissues the analysis was carried out in the same way as for table 5 the possibility that tissues from a particular anatomical site with high dysregulation caused the bias was addressed. oral cavity, tongue, and larynx are the major three sites that together contribute to about 80% of the donated head and neck cancer tissues. cops5-associated go and kegg pathways were determined for the tissues from the particular sites, and confirmed that correlation of cops expression with mitochondria associated genes were not as predominant as in the normal tissues. to further support the finding, enrichment score was calculated for the genes in the oxidative phosphorylation kegg pathway using the rnaseq data of normal head and neck tissues. the obtained es value,, was much higher than the average es value based on randomly assigned gene pools.fig. enrichment score for kegg oxidative phosphorylation pathway. a the enrichment score was calculated for the oxidative phosphorylation kegg pathway. genes are aligned on x-axis based on their correlation coefficients to the cops expression. the hsa pathway contained genes in the entire list. b randomly selected genes were used to calculate es using the r cops gene list as in a. the histogram was generated by iteration of times. the arrow indicates the es for hsa in addition, es values corresponding to all the other kegg pathways were calculated for the matched normal and tumor tissues. the es value of the oxidative phosphorylation pathway found to be among the highest. the ubiquinone biosynthesis pathway generated the highest es value, which is an essential factor in the respiratory chain and requires phenylalanine for its synthesis. thus, the enrichment of genes involved in mitochondria-related pathway based on the cops correlation was exceptionally high.table enrichment scores for the cops gene on kegg pathways all kegg pathways found in the entire gene list of the tcga rna expression database were individually scanned through the r cops alignment, and the maximum as well as minimum were calculated as previously described. the kegg pathways that generated the highest es values are shown. nh: number of genes belonging to a kegg pathway in the entire list. a p value for a particular kegg pathway with an nh was determined using mean and standard deviation obtained by the following method. nk genes were randomly selected from the r cops list and calculated an es value, which was reiterated for times to generate the mean and standard deviation of es to elucidate the significance of the cop coordination with oxidative phosphorylation compared to other genes, we calculated es values relative to the hsa kegg pathway for the entire gene list in the rnaseq data of the matched normal and tumor tissues to calculate fdr. except for the cops7b gene, all the cop subunit genes showed highly significant correlations with the genes in the oxidative phosphorylation pathway. not surprisingly, most of the genes in the hsa oxidative phosphorylation pathway resulted in the high es/low fdr values, validating our approach to detect a linkage of a gene to the oxidative phosphorylation. remarkably, the cop genes were ranked comparably with the hsa genes. the p as well as fdr values were increased for all the cop genes with the tumors, indicating significant decrease in the association of the cop expression to the oxidative phosphorylation related genes. taken all together, we conclude that there is a high synergy between the expression of the cop and mitochondria functions, which implies that the normal physiological function of the cop may be inseparable from the mitochondrial activity. analysis of gene expression of other cancers in the tcga database finally, the question was addressed whether the association of the expressions of cop genes with those of mitochondrial genes was unique in head and neck tissues. to probe this, the rna expression data of tcga for the lung squamous cell carcinoma tissues were analyzed. tumors originating from the lung upper airway epithelia share similarity with hnscc in its etiology and its involvement of smoking as an important risk factor. first, the synergistic expression pattern of cop genes in the normal tissues was examined, and found to be highly synergistic with one another. exceptions were gps and cops7b, which showed lower correlation coefficients to other cop factors. the high positive correlations were again lost in the tumor tissues. these results showed a remarkable coincidence with those of head and neck normal and tumor tissues.table synchronized expression of the cop genes in the matched tissues of normal and lung squamous cell carcinoma the analysis was carried out in the same way as for table  values in bold letters indicate significantly high cosp correlation with other cop components the functional pathway analysis with the entire gene list resulted in the high correlation of the cops expression with the mitochondrial genes. unlike the case of head and neck squamous cell carcinoma, however, the coordinated expression of cop and the mitochondrial genes was maintained well in the lung squamous cell carcinoma tissues. discussion cancer genetics and analyses of the tcga database mainly focus on the loss of genomic integrity. mutations and single nucleotide polymorphism in dna show the high rate of permanent alteration of the genome, epigenetic modification of the genome, rna and mirna expression alterations. tcga data sets are available publicly for these types of analyses. the database is highly reliable and allowed us to examine cancer genomics with more than tumor specimens. in addition to the tumor specimens, data from about to normal tissues are made available for most types of cancers with a notable exception for glioblastoma. although the data sizes are smaller than those for the cancer tissues, they are arguably one of the most comprehensive databases for genomics of human specimens of non-disease status with pathological and clinical scrutiny, which should provide valuable information for elucidating functions of a gene under the normal physiology. we hypothesized that a synergistic expression of genes may be observed in the normal tissues, when expressions of a group of genes need to be coordinated for their proper cellular activities. cop signalosome is an ideal subject to test this possibility, because all eight subunits are essential for the proper cop function; deficiency in one of the cop genes decreased stability of the other subunits and caused lethality. in this study the cop signalosome genes showed highly synergistic expressions in the normal head and neck tissues, except for gps and cops7b. remarkably, these characteristics, including the lower r values for gps and cops7b, were maintained well even in the lung and breast normal tissues. the coordinated expression of the cop genes decreased in the tumor tissues. therefore, this coordination could not have been detected by analyzing tumor tissues. likewise, because all cells cultured in vitro are under non-physiological conditions and thus transformed at some extent, it would have been difficult to obtain convincing evidence for the expression correlation of the cop subunit genes based on in vitro cell biological studies alone. it is interesting to observe the lowest coordination of gps expression compared to the other cop components. besides being a subunit of cop gps is known for its function in suppression of g-protein and mitogen-activated signal transduction in mammalian cells. this cop9-independent function quite likely requires a regulation of the gps gene separate from the other cop genes. cops7b not only had the lowest r value, but also significantly lower expression compared to cops7a. while cops7a and cops7b genes are encoded in different human chromosomes, chromosome and respectively, there is a high homology between cops7a and cops7b in the amino acid sequence. thus, cops7b may have only a minor functional role for cop which may be compensated by cops7a. to our knowledge, using the exact physical map of the transcription start sites available from the ensemble genome database has not been applied previously to elucidate effects of the relative locations of genes on the gene expression correlations. in this study, we used this approach to find the chromosomal locations of the genes as a dominant factor causing the loss of the coordinated expressions of the cop genes. much fewer number of genes were under the control of this cis-effect in case of normal tissues, indicating that most genes in normal tissues are regulated independently of their chromosomal locations and perhaps more based on their functional groups, as supported by this study. therefore, important implications derived from this analysis are that 1) analyzing tumor tissue data may not identify a normal function of a protein family such as cop and 2) the qualities of the normal tissue data available in tcga are high enough to help identifying unknown expressional linkages with functional pathways. although the primary focus of this study was cop expression profiling, we expanded the analysis of the influence of the chromosomal location for all the genes that were analyzable in the rnaseq data, and the data for cancer genes are shown in additional file 7: table s it may be interesting to elucidate whether these proto-oncogenes have tight association with the chromosomal locations in other types of tumor tissues available in tcga. the observation that the expression of the cops gene was synchronized with those involved in mitochondrial function was unexpected. the correlation coefficients of the cops and top ranked genes such as ndufb and mrps were remarkably high with statistical significance. moreover, not only a handful of genes involved in mitochondria but also the total pathway analysis revealed that mitochondria-related cellular functions and pathways were predominantly over-represented in the cops correlation coefficient table. although there have not been reports that show direct function of cop for mitochondria, correlation analyses with the tp of which functions have been extensively studied linked it with known cellular functions including cell cycle and dna damage response pathways. the result validates the approach in this study, and support the prediction that the cop signalosome has a novel functional link to mitochondrial activities. few studies have linked cop to mitochondria-related cellular activities. jab1/cops was previously identified as an apoptosis inducing factor through modulation of mitochondrial bh3-domain containing protein bclgs. a study found that in neurospora crassa, three subunits of the cop were required to recover from deficiency of the alternative oxidase function, which is important for mitochondrial functions. however, no follow up studies have been carried out for mammalian cells. on the other hand, pathways involved in ribosomal functions were associated with the cops expression in the tumor tissues. the significance of this link is currently unclear. however, the cop structure is homologous to translation initiation factor, another high-ordered complex in the cells, and our speculation is that cops may have an overlapping function in the protein synthesis in transformed cells. the quantity and quality of the genomics database will only be improved. the present study underscores the importance and usefulness of examining the genomics of normal tissues. the results would provide novel functional links of a gene, which was not possible by analyzing the data set of pathophysiological specimens. this approach can reasonably be applied to other genes to identify pathways that have not been known due to scarce information for genomics information of normal tissues.  the present study investigated cops gene expression using the tcga database. the analysis revealed the highly synergistic mrna expression pattern among cop genes in normal oral tissues. the coordination was abrogated in the tumor tissues, and the regulatory mechanism was taken over by a cis-acting gene regulation. further analysis revealed an unexpected expression correlation between the cops gene and a number of mitochondria-related genes, postulating the possible functional role of the cop signalosome for mitochondrial activities. the quality of the rna expression data available in tcga is high enough to allow us gene expression analysis in normal tissues to identify gene functions under normal physiology, which may not be recognized by the data obtained from tumor tissues. additional files for analysis of rna expression of normal and cancer tissues reveal high correlation of cop gene expression with respiratory chain complex components by christina a. wicker and tadahide izumi.additional file 1: figure s expressions of cop genes in normal oral tissues. the rna expression levels of cop genes from normal oral tissues are shown. horizontal lines in the boxes: median values. top and bottom whiskers denote the first and the third quartile. outliers and < q - * ) are shown as individual dots. figure s pair-wise plot for expressions of cop genes in the matched lung tissues of normal and squamous cell carcinoma. the plots were generated as described in fig.  except for the data set of the normal and their matched lung squamous cell carcinoma from tcga was used. the number of samples =  additional file 2: table s the let- expression profilesof the matched head and neck normal and tumor tissues. correlation coefficients for individual let- isoforms relative to cop genes are shown. asterisks: p <. table s to support the result of table additional normal head and neck tissues were included to calculate correlation of expressions among the cop genes. correlation coefficients of cop genes versus cops were determined from randomly selected tumor specimens from the entire pool of the calculations were repeated for times to obtain their mean and standard deviation values, which were used to compare the difference from rcops of normal tissues. table s rnaseq data from young and elderly patients were used to calculate the correlation coefficients among the cop genes. table s normal oral cavity tissues were pooled based on the tcga clinical data, and the correlation coefficients and their p values for cops were calculated as described in the main article. table s data based on the normal tongue tissues were pooled based on the clinical data in the tcga database, and the correlation coefficients and their p values for cops were calculated. table s normal tissues from smokers were pooled based on the clinical data in the tcga database, and the correlation coefficients and their p values for cops were calculated. table s nh: number of genes found within mbp of each cop gene; mean, sd, p and q values were determined as described in materials and methods. table s es for cop genes for all kegg pathways. table s es for cop genes for the oxidative phosphorylation pathway. additional file 3: table s the cops correlation coefficient table list includes highest |r| values relative to the expression of cops in the normal and tumor tissues of head and neck. additional file 4: table s go and kegg pathways that are over-represented in the normal and tumor oral tissues based on the cops correlation coefficient table. additional file 5: table s go and kegg pathways over-represented by cops expression correlation in the two major anatomical sites in the head and neck. additional file 6: table s go and kegg pathways over-represented in the normal and tumor lung tissues based on the cops correlation coefficient table. additional file 7: table s genes under the chromosomal cis-effect in the hnscc tissues additional file 8: table s kegg pathway analysis for the tp gene funding this work was supported by the nih training grant t32es and nci grants ca cancer center supporting grant p30ca and ctsa ul1tr availability of data and materials all the data discussed in this study are available as in main figures or in the supplemental data. authors’ contributions ti conceived the project idea and wrote necessary programs. ti and caw performed the data analysis. ti and caw wrote and edited the manuscript, and approved the final manuscript. competing interests the authors declare that they have no competing interests. consent for publication not applicable. ethics approval and consent to participate not applicable.
the results indicate that mrna expression data for the matched normal tissues available in tcga are statistically reliable, and are highly useful to assess novel associations of genes with functional pathways in normal physiology as well as in the cancer tissues. the cop signalosome, composed of eight subunits, is implicated in cancer genetics with its deneddylase activity to modulate cellular concentration of oncogenic proteins such as ikb and tgfβ. pathway analysis revealed a high coordination of the expression of the cops and the other cop genes with mitochondria-related functional pathways, including genes encoding the respiratory chain complex. primarily focusing on gene expression data of the normal tissues of the head and neck, the cancer genome atlas database was used to identify groups of genes that were expressed synergistically with the cop genes, particularly with the cops, which possesses the catalytic activity of cop expressions of seven of the cop genes were found to be highly synergistic in the normal tissues. in contrast, the tumor tissues decreased the coordinated expression pattern of cop genes. this study revealed the significant correlation between the expressions of the cop genes and those related to the mitochondrial activity.
 the cop signalosome, composed of eight subunits, is implicated in cancer genetics with its deneddylase activity to modulate cellular concentration of oncogenic proteins such as ikb and tgfβ. however, its function in the normal cell physiology remains elusive. primarily focusing on gene expression data of the normal tissues of the head and neck, the cancer genome atlas database was used to identify groups of genes that were expressed synergistically with the cop genes, particularly with the cops, which possesses the catalytic activity of cop  expressions of seven of the cop genes were found to be highly synergistic in the normal tissues. in contrast, the tumor tissues decreased the coordinated expression pattern of cop genes. pathway analysis revealed a high coordination of the expression of the cops and the other cop genes with mitochondria-related functional pathways, including genes encoding the respiratory chain complex.  the results indicate that mrna expression data for the matched normal tissues available in tcga are statistically reliable, and are highly useful to assess novel associations of genes with functional pathways in normal physiology as well as in the cancer tissues. this study revealed the significant correlation between the expressions of the cop genes and those related to the mitochondrial activity. electronic supplementary material the online version of this article contains supplementary material, which is available to authorized users. center for clinical and translational scienceul1tr000117izumi tadahide http://dx.doi.org/ /100000054national cancer institute ca098664izumi tadahide http://dx.doi.org/ /100000054national cancer institutep30ca17758izumi tadahide national institute of environmental health sciences t32es007266wicker christina a. issue-copyright-statement© the author 2016 