This paper was converted on www.awesomepapers.org from LaTeX by an anonymous user.
Want to know more? Visit the Converter page.

A Method for Massively Parallel Analysis of Time Series

Yi H. Yan Yi H. Yan: University of Georgia, Institute of Bioinformatics, Athens, 30602, USA Elizabeth D. Trippe Elizabeth D. Trippe: University of Georgia, Department of Mathematics, Athens, 30602, USA, Middle Georgia State University, Department of Mathematics, Macon, 31206, USA  and  Juan B. Gutierrez Juan B. Gutierrez: jgutierr@uga.edu University of Georgia, Department of Mathematics and Institute of Bioinformatics, Athens, 30602, USA
(Date: August 2, 2025; Date: August 2, 2025)
Abstract.

Quantification of system-wide perturbations from time series -omic data (i.e. a large number of variables with multiple measures in time) provides the basis for many down-stream hypothesis generating tools. Here we propose a method, Massively Parallel Analysis of Time Series (MPATS) that can be applied to quantify transcriptome-wide perturbations. The proposed method characterizes each individual time series through its 1\ell_{1} distance to every other time series. Application of MPATS to compare biological conditions produces a ranked list of time series based on their magnitude of differences in their 1\ell_{1} representation, which then can be further interpreted through enrichment analysis. The performance of MPATS was validated through its application to a study of IFNα\alpha dendritic cell responses to viral and bacterial infection. In conjunction with Gene Set Enrichment Analysis (GSEA), MPATS produced consistently identified signature gene sets of anti-bacterial and anti-viral response. Traditional methods such as EDGE and GSEA Time Series (GSEA-TS) failed to identify the relevant signature gene sets. Furthermore, the results of MPATS highlighted the crucial functional difference between STAT1/STAT2 during anti-viral and anti-bacterial response. In our simulation study, MPATS exhibited acceptable performance with small group size (n = 3), when the appropriate effect size is considered. This method can be easily adopted for other -omic data types.

1. Introduction

The advent of high throughput molecular technologies has allowed us to measure a large number of variables simultaneously. Transcriptomic experiments typically measure the abundance of 20,000 or so transcripts. The number of variables measured by other “omic” technologies such as lipidomic, proteomic, glycomic and metabolomic studies are of similar magnitude. Time series -omic data refers to molecular snapshots taken using these -omic technologies on a time trajectory. Time series -omic data captures time dependent molecular dynamics and is a powerful tool in the study of disease progression [2], developmental process [12] and vaccination [11].

A variety of tools have been developed for the analysis of time series -omic data. Methods such as MESS[3] and EDGE[14] aim to discover individual gene expression time series that are significantly different between two experimental conditions. Other methods such as TcGSA [8] ,CAMERA [18] , and GSEA-TS aim to find time series of pre-defined gene sets that are significantly different among groups. Furthermore, clustering based time series analysis tools [15, 5] have also been developed.

Current time series analysis tools focus on finding gene expression time series or pre-defined gene sets that most likely have changed between groups, but ignore the changes in pairwise gene dynamics. Pair wise dynamics between genes can be quantified using correlation, mutual information or a distance metric. Differential correlation has been used to study gene association with the clinical outcome of lung cancer [13] and estrogen receptor modulation in hormonal cancers [10]. In this paper, we propose a novel method, Massive Parallel Analysis of Time Series (MPATS), for the analysis of time series data based on the detection of differential pairwise 1\ell_{1} distances. MPATS is based on the use of a linear mixed model for gene expression time series. By characterizing each gene expression time series through its 1\ell_{1} distance to all other time series within the system, MPATS allows for the clear quantification of the impact of perturbation on each time series in the context of the biological system.

A simulation study was conducted to demonstrate the statistical performance of our method and a motivational study was used to validate the biological relevance of the analysis result produced by MPATS. In the motivational study, time series transcriptomic data were used to characterize the responses of IFNα\alpha dendritic cells to different antigens [1]. In contrast with traditional time series analysis methods such as EDGE and GSEA-TS, MPATS was able to identify signature gene modules that distinguish between anti-viral and anti-bacterial responses of IFNα\alpha dendritic cells.

This paper is organized as follows: section 2 contains the derivations, theoretical justifications and description of our method. Section 3 provides the assessment of statistical power of our method based on simulated data. Section 4 contains the application of our method to a recent study of IFNα\alpha dendritic cells. Section 5 provides discussion of the performance of our method and future directions.

2. Method Description

2.1. Linear Mixed Model of Time Series

Each individual time series is denoted Gij(t)G_{ij}(t) and represents the abundance of gene ii from individual jj at time point tt, where i{1,2,n}i\in\{1,2,...n\}, j{1,2,m}j\in\{1,2,...m\} and t{1,2h}t\in\{1,2...h\}. We used a mixed linear model [14] for the time series where:

Gij(t)=gij(t)+νt+ϵ.G_{ij}(t)=g_{ij}(t)+\nu_{t}+\epsilon. (1)

The above model implies that the observed entity Gij(t)G_{ij}(t) can be explained by the additive effect of its time dependent mean response gij(t)g_{ij}(t), individual and time based variation νt\nu_{t}, and instrumentation error ϵ\epsilon. We assume that both νt\nu_{t} and ϵ\epsilon follow Gaussian distributions, νt𝒩(0,νt)\nu_{t}\sim\mathcal{N}(0,\nu_{t}) and ϵ𝒩(0,ϵ)\epsilon\sim\mathcal{N}(0,\epsilon). Equation (1) can be further simplified to:

Gij(t)𝒩(gij(t),νt+ϵ).G_{ij}(t)\sim\mathcal{N}\left(g_{ij}(t),\nu_{t}+\epsilon\right). (2)

2.2. 1\ell_{1} Distance Between Two Time Series

Let 1(p,q,j)\ell_{1}(p,q,j) denote the 1\ell_{1} distance between the mean response curves of genes pp and qq for subject jj.

1(p,q,j)=t=1h|gpj(t)gqj(t)|.\ell_{1}(p,q,j)=\sum_{t=1}^{h}|g_{pj}(t)-g_{qj}(t)|. (3)

Let X(p,q,j)X(p,q,j) be the observed 1\ell_{1} distance between the two genes pp and qq for subject jj.

X(p,q,j)=t=1hx(p,q,j,t),X(p,q,j)=\sum_{t=1}^{h}x(p,q,j,t), (4)

where

x(p,q,j,t)=|Gpj(t)Gqj(t)||𝒩(μt,σt2)|,x(p,q,j,t)=|G_{pj}(t)-G_{qj}(t)|\sim|\mathcal{N}\left(\mu_{t},\sigma_{t}^{2}\right)|, (5)

and

μt=gpj(t)gqj(t),σt2=2(νt+ϵ)\begin{split}\mu_{t}&=g_{pj}(t)-g_{qj}(t),\\ \sigma_{t}^{2}&=2(\nu_{t}+\epsilon)\end{split} (6)

Through numerical investigation, Tsagris et al [17] saw that for any folded normal distribution FF, where

F|𝒩(μ,σ2)|,F\sim|\mathcal{N}(\mu,\sigma^{2})|, (7)

if μ3σ2\mu\geq 3\sigma^{2}, then FF can be well approximated with a normal distribution with mean |μ||\mu| and variance σ2\sigma^{2}.

Assuming that for any tt, |μt|3σt2|\mu_{t}|\geq 3\sigma_{t}^{2}, then

X(p,q,j)t=1h𝒩(μt,σt2).X(p,q,j)\sim\sum_{t=1}^{h}\mathcal{N}(\mu_{t},\sigma_{t}^{2}). (8)

Which means

X(p,q,j)𝒩(t=1h|μt|,t=1hσt2).X(p,q,j)\sim\mathcal{N}(\sum_{t=1}^{h}|\mu_{t}|,\sum_{t=1}^{h}\sigma_{t}^{2}). (9)

According to equation (6) and (3) t=1h|μt|=1(p,q,j)\sum_{t=1}^{h}|\mu_{t}|=\ell_{1}(p,q,j), and

X(p,q,j)𝒩(1(p,q,j),t=1hσt2).X(p,q,j)\sim\mathcal{N}(\ell_{1}(p,q,j),\sum_{t=1}^{h}\sigma_{t}^{2}). (10)

Which implies that given sufficiently small variance for all tt, the observed 1\ell_{1} distance between two gene expression time series follows an approximate normal distribution with mean equal to the true 1\ell_{1} distance.

In the case that |μt|<3σt2|\mu_{t}|<3\sigma_{t}^{2}, the mean of X(p,q,j)X(p,q,j) does not necessarily equal to 1(p,q,j)\ell_{1}(p,q,j).

Let cc be a constant, such that

c>max10σt2t{1,2,3h}.c>\max 10\sigma_{t}^{2}\quad\forall\quad t\in\{1,2,3...h\}. (11)

Then let

y1(p,q,t,c,j)=Gpj(t)Gqj(t)+c,|𝒩(μ1,σt2)|,\begin{split}y_{1}(p,q,t,c,j)&=G_{pj}(t)-G_{qj}(t)+c,\\ &\sim|\mathcal{N}\left(\mu_{1},\sigma_{t}^{2}\right)|,\end{split} (12)

where

μ1=gpj(t)gqj(t)+c.\mu_{1}=g_{pj}(t)-g_{qj}(t)+c. (13)

Let

y2(p,q,t,c,j)=Gqj(t)Gpj(t)+c,|𝒩(μ2,σt2)|,\begin{split}y_{2}(p,q,t,c,j)&=G_{qj}(t)-G_{pj}(t)+c,\\ &\sim|\mathcal{N}\left(\mu_{2},\sigma_{t}^{2}\right)|,\end{split} (14)

where

μ2=gqj(t)gpj(t)+c.\mu_{2}=g_{qj}(t)-g_{pj}(t)+c. (15)

Because c>max10σt2t{1,2,3h}c>\max 10\sigma_{t}^{2}\quad\forall\quad t\in\{1,2,3...h\},

y1(p,q,t,c,j)𝒩(μ1,σt2),y2(p,q,t,c,j)𝒩(μ2,σt2).\begin{split}y_{1}(p,q,t,c,j)&\sim\mathcal{N}(\mu_{1},\sigma_{t}^{2}),\\ y_{2}(p,q,t,c,j)&\sim\mathcal{N}(\mu_{2},\sigma_{t}^{2}).\\ \end{split} (16)

Furthermore,

max(μ1,μ2)=|gpj(t)gqj(t)|+c.\max(\mu_{1},\mu_{2})=|g_{pj}(t)-g_{qj}(t)|+c. (17)

Let

z(p,q,t,c,j)={y1(p,q,t,c,j)if μ1μ2,y2(p,q,t,c,j)if μ1<μ2.z(p,q,t,c,j)=\begin{cases}y_{1}(p,q,t,c,j)&\mbox{if }\mu_{1}\geq\mu_{2},\\ y_{2}(p,q,t,c,j)&\mbox{if }\mu_{1}<\mu_{2}.\end{cases} (18)

Let

Z(p,q,c,j)=t=1hz(p,q,t,c),Z(p,q,c,j)=\sum_{t=1}^{h}z(p,q,t,c), (19)

then,

Z(p,q,c,j)𝒩(1(p,q,j)+hc,t=1hσt2).Z(p,q,c,j)\sim\mathcal{N}(\ell_{1}(p,q,j)+hc,\sum_{t=1}^{h}\sigma_{t}^{2}). (20)

Which implies that Z(p,q,c,j)Z(p,q,c,j) follows a normal distribution with mean equal to the true 1\ell_{1} distance plus some known constan t.

2.3. Detection of Significantly Changed pq\ell_{pq} Between Groups

Each individual gene expression time series can be characterized by its 1\ell_{1} distances to every other gene expression time series. A perturbation of the system (e.g. disease progression) shifts the mean response curve of each gene expression time series with an unknown magnitude and direction.

The impact of the perturbation on the system can be quantified by the number of 1\ell_{1} distances that have changed due to perturbation. The impact of the perturbation on each individual gene time series can be quantified by the number of pairwise distances to this gene that have changed.

In the case of a two-group experiment, we define two sets of subjects, where S1={1,2,3,4a}S_{1}=\{1,2,3,4...a\} and S2={a+1,a+2m}S_{2}=\{a+1,a+2...m\}. For a given pair of entities pp and qq, let Apqc={Z(p,q,j,c):jS1}A_{pqc}=\{Z(p,q,j,c):j\in S_{1}\} and Bpqc={Z(p,q,j,c):jS2}B_{pqc}=\{Z(p,q,j,c):j\in S_{2}\}. As shown in (20), both ApqcA_{pqc} and BpqcB_{pqc} follow a normal distribution with mean corresponding to the 1\ell_{1} distance between gene pp and qq within each group. Welch’s test is used to determine if ApqcA_{pqc} and BpqcB_{pqc} have the same mean.

Let P(p,q)P(p,q) be equal to the p-value from Welch’s test. The results of conducting Welch’s test for every pair of pp and qq can be captured in a nn by nn matrix MM, where

Mpq={1,P(pq)<δ0,P(pq)δ,M_{pq}=\begin{cases}1,P(pq)<\delta\\ 0,P(pq)\geq\delta,\end{cases}

where δ\delta is the threshold p value.

For a gene pp, its contribution to the perturbation can be quantified by its perturbation score (PS)(PS), where:

PS=q=1nMpq.PS=\sum_{q=1}^{n}M_{pq}. (21)

The genes are then ranked according to their perturbation score from highest to lowest.

3. Power Assessment Using Simulated Data

A simulation study was conducted to evaluate the statistical power of MPATS. The simulation framework was chosen based on the motivational study. For each simulation setup, 5000 gene expression time series were simulated and each time series contained 5 time points. The simulated values of each gene at each time point were drawn from a normal distribution with mean and variance estimated from the original data. 10 percent of the gene expression time series were perturbed in the treatment group. Receiver operating characteristic (ROC) curves were generated to explore the statistical performance of the proposed method for different group sizes and effect sizes. The effect size (ES) refers to the true difference in 1\ell_{1} distance, and the category of effect sizes can be found in Table 1. Statistical performance of the proposed method are presented in Figure 1. At a fixed specificity of 90%90\%, MPATS has statistical power greater than 0.70.7 with n=3n=3 for medium and high effect sizes. With n=10n=10, MPATS can reach a power greater than 0.70.7 with a fixed specificity of 90%90\% for low effect size.

Effect Size Description

Category ESES
Minimal 0<ES0.250<ES\leq 0.25
Low 0.25<ES0.50.25<ES\leq 0.5
Medium 0.5<ES10.5<ES\leq 1
High 1<ES1<ES
Table 1. Categories of Effect Size

ROC Curve of Different Sample Size and Effect Size

Refer to caption

Figure 1. ROC curves were generated for 4 different effect sizes. For each effect size, three different sample sizes were tested, n=3n=3, n=5n=5 and n=10n=10. For each simulation, 5000 gene expression time series were generated.

4. Experimental Results

4.0.1. Motivating Study: IFNα\alpha dendritic cell response to antigen

As part of a study to understand the different responses of dendritic cells to various antigens, microarray time series data were collected for IFNα\alpha dendritic cells and IL4 dendritic cells challenged with different antigens. Data were collected at (1,2,6,12,24) hours after challenge. At each time point, three biological replicates were taken. The detailed description of the study can be found in Banchereau et al [1].

For the purpose of the present paper, we focused on the time series of IFNα\alpha dendritic cells grown on media alone (MEDIA), challenged with H1N1 virus (H1N1) and heat killed Staphylococcus aureus (HKSA). The data used in this study can be found on Gene Omnibus (GSE44720).

4.0.2. Application of MPATS, EDGE, and GSEA-Time Series

MPATS and EDGE were used to conduct pairwise comparisons among the two treatment groups and the control group. GSEA-TS was applied to the time series data of the two treatment groups alone. During the application of MPATS, a q-value of 0.250.25 was used as the threshold to determine whether a pairwise 1\ell_{1} distance had changed significantly between groups. The ranked gene list produced by MPATS was analyzed using GSEA preranked [16]. The results of EDGE consisted of q-values for each gene. A ranked list was produced by ranking the genes according to their q-values in ascending order. The ranked gene list produced by EDGE was also analyzed using GSEA preranked. The gene sets tested for enrichment included: GSEA Hallmark gene sets, GO biological process gene sets, KEGG gene sets, and the signature gene modules described by Banchereau et al [1].

4.1. Clustering of pairwise 1\ell_{1} distance representation of samples

Hierarchical clustering of the samples based on their pairwise 1\ell_{1} distance representation reveals differences in biological conditions (Fig 2). Furthermore, projection of the data into the first and second principle component space also demonstrate clustering based on biological condition (Fig 2). The hierarchical clustering based on biological condition can be reproduced using only the top 100 1\ell_{1} distances with the largest variance across all three biological conditions (Fig 3).

Refer to caption
Figure 2. Both hierarchical clustering of pairwise 1\ell_{1} distance and projection into principle component space reveal underlying biological conditions.
Refer to caption
Figure 3. Top 100 1\ell_{1} distances with largest variance can be used to cluster the samples based on biological condition. Each row refers to the 1\ell_{1} distance between two genes. Each column represents a sample.
Refer to caption
Figure 4. Percentage of genes associated with at least one 1\ell_{1} distance smaller than a threshold (q).

Only a subset of genes are associated with 1\ell_{1} distances that have a small q-value. For the comparison between the H1N1 group and the MEDIA group, only 10%10\% of the genes are associated with a 1\ell_{1} distance with a q-value smaller than 0.0750.075 (Figure 4). Furthermore, the highly significantly changed 1\ell_{1} distances are associated with a small portion of the genes as shown in Figure 5. For the comparison the between the H1N1 group and the MEDIA group, 70%70\% of the 1\ell_{1} distances with qvalue0.1q-value\leq 0.1 are associated with 20%20\% of the genes. This non-uniform distribution of these significantly changed 1\ell_{1} distances provides the basis for downstream enrichment analysis. An example of a significantly changed 1\ell_{1} distance between two genes is shown in Figure 6. The ranking of each gene is based on the number of significantly changed 1\ell_{1} distance associated with that gene. Figure 7 shows that he ranking of each gene is strongly associated with the change in the time series of the gene itself.

Refer to caption
Figure 5. Highly significantly changed 1\ell_{1} distance are associated with a small sub portion of the total genes.
Refer to caption
Figure 6. 1\ell_{1} distance between ACAT1 and LYSE changed significantly (q0.005)(q\leq 0.005) between the H1N1 challenged group and HKSA challenged group. Time series from all three biological replicated are displayed.
Refer to caption
Figure 7. Examples of most significantly changed 1\ell_{1} distances associated with genes of different ranking. The red solid lines are the time series expression profile of the gene of interest in each row. LY6E has rank 1, OAS3 has rank 100 and CCDC90A has rank 1000. The ranking of genes strongly associates with the changes in temporal pattern.

4.1.1. MPATS Identifies Signature Gene Sets

Number of Perturbed Gene Sets Discovered by Each Method

HKSA H1N1 H1N1\SA
MPATS 326 348 270
EDGE 24 39 26
GSEA-TS 48 42 NA
Table 2. Three pairwise comparisons were conducted. Comparison between control group and treatment group challenged by SA (SA). Comparison between control group and treatment group challenged by H1N1 (H1N1). Comparison between the two treatment groups (H1N1\SA).

MPATS identified the greatest number of perturbed gene sets out of all three methods. The number of perturbed gene sets identified by each method are reported in Table 2.The differences in the number of perturbed gene sets identified among the three methods highlight the sensitivity of MPATS. This increased sensitivity is due to the fact that MPATS ranks the genes based on their contribution to system-wide dynamic differences between biological conditions. The large number of gene sets perturbed highlights the tremendous impact of antigen challenge on the cell culture. A majority of gene sets and pathways identified are interferon response related pathways. The enrichment results are included in the supplementary material.

In the original study, a list of signature modules were identified for INFα\alpha cell response to different antigens. Enrichment analysis of both MPATS and EDGE results captured 9 out of 10 signature modules of INFα\alpha response to H1N1 challenge, where GSEA-TS only captured 5. All three method performed worse for the analysis of signature response modules of anti-bacterial response, capturing only 30%30\% of the signature modules. This could be attributed to the low signal strength of these modules themselves in comparison to the signal strength of the H1N1 modules.

Number of Enriched Signature Modules

SA H1N1
MPATS 11/30 9/10
EDGE 9/30 9/10
GSEA-TS 12/30 5/10
Table 3. MPATS results are enriched in the signature modules of INFα\alpha cells challenged with different antigens.

The analysis comparing IFNα\alpha response to H1N1 and HKSA are enriched in interferon related gene sets, such as response to interferon gamma, TNF signal of beta kappa B and MTORC signaling. Analysis of the top 30 genes using String10 shows a enrichment in protein protein interaction (PPI) with pvlaue<0.00001p-vlaue<0.00001. Furthermore, the PPIs are clustered around STAT1 and STAT2 (Figure 8). This result agrees with recent findings of the crucial differences in the role of STAT1/STAT2 in anti-viral and anti-bacterial responses [4, 7, 6, 9].

Refer to caption
Figure 8. PPI Map of the top 30 genes for the comparison of H1N1 and HKSA.

5. Discussion

In this paper, we presented a novel framework to quantify the perturbation of time series for two group experiments. The quantification of perturbation of time series provides context for downstream functional analysis. We were able to apply this framework to explore anti-viral and anti-bacterial response of IFNα\alpha dendritic cells. We found specific pathways that differentiate between viral and bacterial challenge of IFNα\alpha dendritic cells and identified STAT1/STAT2 as important regulatory elements differentiating these responses.

MPATS quantifies system wide perturbations and individual gene perturbations by characterizing each time series by its 1\ell_{1} distance to every other time series. The biological significance of this intuitive method was demonstrated through analysis of the motivational study. The top ranked genes are not only enriched in PPIs but are also enriched in signature gene sets. The performance of MPATS remained stable through all three analyses, whereas EDGE and GSEA-TS identified varying numbers of perturbed gene sets. This is probably due to the fact that ranking by p-value of change is not biologically informative because a change of small magnitude with low variance can generate a low p-value. Ranking time series using variance normalized effect size should produce results similar to MPATS. The framework of MPATS is based on a linear mixed model and assumes normally distributed variance, and the use of pairwise dynamics to quantify magnitude of change for each individual entity can be easily expanded to other -omic data types.

In conclusion, MPATS complements existing time series analysis methods by providing a more biologically informative ranked list of genes of interest in addition to detecting genes that have experienced large perturbations but are overlooked by existing methods.

Acknowledgments

Research reported in this publication was supported by a National Institute of Allergy and Infectious Diseases cooperative agreement of the National Institutes of Health under award number HHSN272201200031C.

References

  • [1] R. Banchereau, N. Baldwin, A.-M. Cepika, S. Athale, Y. Xue, I. Y. Chun, P. Metang, A. Cheruku, I. Berthier, I. Gayet, et al. Transcriptional specialization of human dendritic cell subsets in response to microbial vaccines. Nature communications, 5, 2014.
  • [2] C. Bécavin, N. Tchitchek, C. Mintsa-Eya, A. Lesne, and A. Benecke. Improving the efficiency of multidimensional scaling in the analysis of high-dimensional data using singular value decomposition. Bioinformatics, 27(10):1413–1421, 2011.
  • [3] M. Berk, C. Hemingway, M. Levin, and G. Montana. Longitudinal analysis of gene expression profiles using functional mixed-effects models. In Advanced Statistical Methods for the Analysis of Large Data-Sets, pages 57–67. Springer, 2012.
  • [4] K. Blaszczyk, A. Olejnik, H. Nowicka, L. Ozgyin, Y.-L. Chen, S. Chmielewski, K. Kostyrko, J. Wesoly, B. L. Balint, C.-K. Lee, et al. Stat2/irf9 directs a prolonged isgf3-like transcriptional response and antiviral activity in the absence of stat1. Biochemical Journal, 466(3):511–524, 2015.
  • [5] J. Ernst, G. J. Nau, and Z. Bar-Joseph. Clustering short time series gene expression data. Bioinformatics, 21(suppl 1):i159–i168, 2005.
  • [6] B. Hahm, M. J. Trifilo, E. I. Zuniga, and M. B. Oldstone. Viruses evade the immune system through type i interferon-mediated stat2-dependent, but stat1-independent, signaling. Immunity, 22(2):247–257, 2005.
  • [7] S. Hambleton, S. Goodbourn, D. F. Young, P. Dickinson, S. M. Mohamad, M. Valappil, N. McGovern, A. J. Cant, S. J. Hackett, P. Ghazal, et al. Stat2 deficiency and susceptibility to viral illness in humans. Proceedings of the National Academy of Sciences, 110(8):3053–3058, 2013.
  • [8] B. P. Hejblum, J. Skinner, and R. Thiébaut. Time-course gene set analysis for longitudinal gene expression data. PLoS Comput Biol, 11(6):e1004310, 2015.
  • [9] M. J. Hofer, W. Li, P. Manders, R. Terry, S. L. Lim, N. J. King, and I. L. Campbell. Mice deficient in stat1 but not stat2 or irf9 develop a lethal cd4+ t-cell-mediated disease following infection with lymphocytic choriomeningitis virus. Journal of virology, 86(12):6932–6946, 2012.
  • [10] T.-H. Hsiao, Y.-C. Chiu, P.-Y. Hsu, T.-P. Lu, L.-C. Lai, M.-H. Tsai, T. H.-M. Huang, E. Y. Chuang, and Y. Chen. Differential network analysis reveals the genome-wide landscape of estrogen receptor modulation in hormonal cancers. Scientific reports, 6, 2016.
  • [11] R. E. Palermo, L. J. Patterson, L. D. Aicher, M. J. Korth, M. Robert-Guroff, and M. G. Katze. Genomic analysis reveals pre-and postchallenge differences in a rhesus macaque aids vaccine trial: insights into mechanisms of vaccine efficacy. Journal of virology, 85(2):1099–1116, 2011.
  • [12] S. Roy, J. Ernst, P. V. Kharchenko, P. Kheradpour, N. Negre, M. L. Eaton, J. M. Landolin, C. A. Bristow, L. Ma, M. F. Lin, et al. Identification of functional elements and regulatory circuits by drosophila modencode. Science, 330(6012):1787–1797, 2010.
  • [13] K. Shedden and J. Taylor. Differential correlation detects complex associations between gene expression and clinical outcomes in lung adenocarcinomas. In Methods of Microarray Data Analysis, pages 121–131. Springer, 2005.
  • [14] J. D. Storey, W. Xiao, J. T. Leek, R. G. Tompkins, and R. W. Davis. Significance analysis of time course microarray experiments. Proceedings of the National Academy of Sciences of the United States of America, 102(36):12837–12842, 2005.
  • [15] J. Straube, A.-D. Gorse, B. E. Huang, K.-A. Lê Cao, et al. A linear mixed model spline framework for analysing time course ‘omics’ data. PloS one, 10(8):e0134540, 2015.
  • [16] A. Subramanian, P. Tamayo, V. K. Mootha, S. Mukherjee, B. L. Ebert, M. A. Gillette, A. Paulovich, S. L. Pomeroy, T. R. Golub, E. S. Lander, et al. Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles. Proceedings of the National Academy of Sciences, 102(43):15545–15550, 2005.
  • [17] M. Tsagris, C. Beneki, and H. Hassani. On the folded normal distribution. Mathematics, 2(1):12–28, 2014.
  • [18] D. Wu and G. K. Smyth. Camera: a competitive gene set test accounting for inter-gene correlation. Nucleic acids research, 40(17):e133–e133, 2012.