Clinical profile and gene expression patterns for mild and severe COVID-19 patients during 28 days. A. Longitudinal sampling schedule for severe and mild COVID-19 donors of peripheral blood (22 samples in total from 8 donors). Three sampling times (black dots) are displayed with respect to the recruitment day (D0). In addition, panel A shows the diagnosis day with positive PCR (green star), symptoms onset day (orange star), and hospitalization day (red cross). B. The COVID-19 progression according to the WHO ordinal scale describes the temporal disease severity for each donor. Severe patients (S1-S4) are displayed in lines with colors scaling from green to red, while mild patients (M1-M4) are shown in blue line colors. The X-axis represents the days relative to the recruitment day. The vertical dot lines indicate the three sampling times used in this study (days 0, 7, and 28). C. Principal Component Analysis plot based on gene expression profiles for mild (blue) and severe (red) COVID-19 patients and grouped by their sampling time (0, 7, and 28 days after recruitment). In addition, peripheral blood samples from two healthy donors are shown in green. PC = Principal Component. D. Heatmap of the 100 most significant differentially expressed genes related to the COVID-19 disease progression. At the bottom, each column corresponds to the sampling points (0, 7, and 28 days since recruitment) of mild and severe patients. Genes are displayed as horizontal rows and are clustered by the similarity of expression profiles, represented by the dendrogram to the left of the heatmap. Red indicates higher expression, while blue means lower expression represented by the z-score of normalized read counts. Some COVID-19 severity-associated genes previously reported are indicated to the right of the heatmap in red.

Clinical characteristics of the cohort

Gene Set Enrichment Analysis (GSEA) shows biological processes and enriched pathways associated to innate and adaptive immune response after SARS-CoV-2 infection. (A – B) Bubble plots show the biological processes (BP) of enriched genes for DEGs between mild and severe for D0 and D7 after recruitment, respectively. (C – D). Bubble plots show the KEGG enriched pathways of enriched genes for DEGs between mild and severe for D0 and D7, respectively. The Count (black circles) represents the number of genes included on each set. Generatio is the ratio between numbers of genes found in the set and total genes of set. The scale-color bar indicates the p-adjustment of each BP or KEGG pathway.

Biological processes and KEGG pathway for genes with differential expression levels over time in mild versus severe patients. Bar plots show Gene ontology analysis for biological processes in light-blue/dark-blue and KEGG pathway in light-green/dark-green, using upregulated genes (A - C - E) and downregulated genes (B – D - F). The size of each bar is according to its -log10 p-value, and name pathways with asterisk indicates a q-value ≤ 0.05.

Gene function network for Natural Killer cell hub-genes with differential expression levels between mild and severe patients during acute phase of COVID-19. (A) KEGG pathway of Natural Killer cell mediated cytotoxicity represents the set NK cell hub-genes upregulated (red-boxes) in mild versus severe patients. The green boxes correspond to genes without differential expression. (B) Heatmap shows the differential expression levels of NK cell hub-genes over time (D0, D7, and D28 after recruitment) separated by mild and severe groups. The expression levels are represented by the z-score of normalized counts. Dendrogram shows the hierarchical clustering of genes. (C) Protein-protein interaction (PPI) network for upregulated genes during the acute phase in mild patients. The network corresponds to the principal clusters with more interaction between proteins and highlight the three most represented pathways: Cytokine-cytokine receptor interaction (blue); NK cell mediated cytotoxicity (red); and Th1 and Th2 cell differentiation (yellow). (D) Time course expression levels for the main protein nodes identified in PPI network during the acute phase of COVID-19. The trajectories of these genes are graphed as days after recruitment (D0, D7, and D28) for mild (red triangle) and severe (green circle) groups and their enrichment is represented by the z-score of normalized counts.

Gene co-expression network analysis among the longitudinal transcriptomic profiling. (A) Gene hierarchical clustering dendrogram of ten detected modules based on Topological Overlap Matrices (TOM) measure. The branches and color bands represent the assigned module. (B) Module-trait relationships (MTRs) between detected modules and disease severity of COVID-19. MTRs are obtained by calculating Pearson correlation between the traits and the module eigengenes. The red and blue colors indicate strong positive or negative correlations, respectively. Rows represent module eigengene (ME) and columns indicate the disease severity of COVID-19, where correlations from each patient are shown from left to right according to their sampling time (D0, D7, and D28). (C) GO enrichment analysis of genes in the blue, brown, and turquoise modules. The color in the bar graphs refers to the module eigengene (ME). Enrichment results are sorted by - log10(p-value) (higher on top) of each biological process GO term. (D-E-F) Network visualization of enriched pathways based on co-expression genes from blue, brown, and turquoise modules, respectively.