Atopic Dermatitis Transcriptome Under Dupilumab and Cyclosporine

Treatment Convergence Hierarchy

Dupilumab vs cyclosporine correlation at each organizational level. The monotonic increase from gene (r = 0.56) to module (r = 0.95) reveals that different gene targets converge into shared biological programs at higher organizational levels.

Novel & Extended Findings

Sample Distribution

Sample counts by skin type, timepoint, and therapy arm. Baseline-only patients (gray) have no therapy annotation since they were not followed up. Three patients contributed only lesional (AL) biopsies at baseline.

Cohort Summary

Library Size Distribution

Samples ordered by total library size (read count). The outlier sample (Patient_31_AL_m0) was removed from downstream analyses. Library sizes range from 11.3M to 43.1M reads (median 22.9M).

UMAP Projection

UMAP of 165 samples (post-QC) colored by skin type. Lesional (AL, red) and non-lesional (AN, blue) samples show partial separation, consistent with the modest silhouette score for skin_type (0.072). Hover for sample details including patient, timepoint, and therapy.

Volcano Plot: AL vs AN at Baseline

Each point represents one of 26,843 tested genes. Red = upregulated in lesional skin (483 genes), blue = downregulated (273), gray = not significant. Top 20 genes labeled in gold. Dashed lines mark thresholds: FDR < 0.05 and |log₂FC| > 1. Model: LMM with patient as random effect (n = 106 samples from 53 patients).

DEG Counts by Comparison

Stacked bars show upregulated (red) and downregulated (blue) gene counts for each differential expression comparison. Cyclosporine yields paradoxically more DEGs (794 AL) than dupilumab (416 AL) despite 2.6× fewer patients, reflecting broader but potentially less reliable transcriptomic impact with only 8 patients.

Gene Set Overlap (Jaccard Similarity)

6×6 Jaccard similarity heatmap between all DE comparisons. Low overlap between dupilumab and cyclosporine AL responses (J = 0.06) contrasts with the moderate overlap between baseline disease signature and dupilumab response (J = 0.11). Sequential blue color scale; darker = less overlap.

Treatment Convergence at Gene Level

Despite only 6% gene overlap (Jaccard = 0.06), the effect direction correlation of r = 0.56 in lesional skin indicates convergent biology — both treatments tend to move the same genes in the same direction, but rarely reach statistical significance for the same targets. Cyclosporine's larger DEG count (794 vs 416) reflects broader but shallower transcriptomic impact, likely amplified by its small sample size (n=8) reducing FDR stringency. The much lower correlation in non-lesional skin (r = 0.28) suggests treatments diverge more where less shared AD signal exists. Dupilumab reverses 20.5% of baseline AL-upregulated genes vs 8.3% for cyclosporine, indicating dupilumab is more targeted at the core AD signature while cyclosporine acts through alternative pathways.

Only 15 genes (3.1% of the 483 baseline AL-upregulated) are reversed by both dupilumab and cyclosporine, defining the minimal shared therapeutic core: S100A12, S100A7A, CCL20, IL36A, SPRR2A/B/D/F, APOBEC3A, ANGPTL4, ADGRF1, CDH3, CLCN1, HEPHL1, VNN3P. These encompass innate immune alarmins (S100 family), chemokines (CCL20), interleukin signaling (IL36A), and epidermal hyperproliferation markers (SPRR family) — representing the irreducible therapeutic target regardless of mechanism. Evidence from Disease Signature and Treatment Response analyses (T2_S1–T2_S5).

Hallmark Pathways: Baseline AL vs AN

All 50 Hallmark pathways ranked by NES. Positive NES (red) = enriched in lesional skin; negative (blue) = enriched in non-lesional. Faded bars are not significant (FDR ≥ 0.25). Top enriched: IFN-gamma response (+2.67), allograft rejection (+2.53), IFN-alpha (+2.43). 37/50 significant at FDR < 0.25.

Treatment Comparison: Pathway-Level Effects

NES values for baseline, dupilumab, and cyclosporine across all 50 Hallmark pathways. Diverging red–blue scale: red = enriched in AL (baseline) or upregulated at m3 (treatment), blue = reverse. Both treatments reverse the baseline pattern, with 27 concordant and 0 discordant pathways.

Per-Sample Pathway Activity

AN_m0 (54 samples) · AN_m3 (26 samples) · AL_m0 (56 samples) · AL_m3 (29 samples)

Individual patient pathway scores (ULM method) for 165 samples × 10 selected pathways. Diverging color scale (red = high, blue = low). Use the vertical slider to scroll through samples. Per-sample treatment convergence r = 0.903, exceeding the GSEA NES correlation.

Pathway Convergence

Despite only 6% gene overlap (Jaccard = 0.06), pathway-level convergence between dupilumab and cyclosporine reaches r = 0.84 (GSEA NES), with all 27 shared significant pathways concordant and zero discordant. Per-sample pathway scoring (ULM) pushes convergence even higher to r = 0.90, adding individual patient resolution. Dupilumab reverses 69.7% of baseline-enriched pathways while cyclosporine reverses 75.8%, both far exceeding their gene-level reversal rates (20.5% and 8.3% respectively). This demonstrates that different genes within shared pathways produce functionally equivalent outcomes.

Cyclosporine's gene-to-pathway reversal amplification is 9.1x (8.3% gene → 75.8% pathway) compared to dupilumab's 3.4x (20.5% → 69.7%). Cyclosporine achieves similar pathway normalization through many coordinated sub-threshold gene changes — consistent with calcineurin inhibition's broader mechanism affecting multiple T-cell lineages simultaneously, producing widespread but individually modest gene perturbations that aggregate into coherent pathway effects. Evidence from Pathway Enrichment analysis (T3_S1–T3_S3).

Top 30 Differential Transcription Factors (AL vs AN)

Top 30 TFs ranked by activity difference (AL − AN). All are AL-enriched (red), reflecting the inflammatory TF landscape of lesional skin. RELA (+1.73) and NFKB1 (+1.62) dominate, confirming NF-kB as the master regulator. STAT3 (+1.31), STAT1 (+1.16), and STAT6 (+1.02) form the JAK-STAT signaling triad.

Treatment TF Activity Changes: Dupilumab vs Cyclosporine

Each point represents one of 761 TFs. x-axis = dupilumab-induced activity change; y-axis = cyclosporine-induced change. Key TFs labeled in gold. The overall correlation (r = 0.768) reflects convergent TF suppression despite distinct mechanisms. Note the subtle divergence: STAT6 moves more with dupilumab, STAT3 more with cyclosporine.

Showing top 30 of 492 significant TFs (FDR < 0.05)

Master Regulators of the AD Transcriptome

NF-kB (RELA, NFKB1) is the dominant upstream regulator of AD lesional skin, driving cytokines, chemokines, and adhesion molecules. The STAT triad — STAT1 (IFN-gamma), STAT3 (IL-6/IL-23/Th17), and STAT6 (IL-4/IL-13/Th2) — reflects the mixed immune activation in chronic AD. Both treatments suppress NF-kB at similar magnitudes (dupilumab ΔRELA = -0.89; cyclosporine ΔRELA = -0.92), establishing convergent TF-level therapeutic suppression. The AP-1 complex (JUN, FOS) and myeloid regulator SPI1/PU.1 confirm the inflammatory and immune cell infiltration signature of lesional skin. TF-level treatment correlation (r = 0.77) sits between gene-level (r = 0.56) and pathway-level (r = 0.84) in the convergence hierarchy.

Dupilumab preferentially suppresses STAT6 (-0.68) over STAT3 (-0.61), while cyclosporine preferentially suppresses STAT3 (-0.78) over STAT6 (-0.58). This provides a TF-level mechanistic explanation: dupilumab blocks IL-4Rα → STAT6 (type 2 immunity), while cyclosporine inhibits calcineurin → NFAT/STAT3 (IL-6/Th17 pathway). The preferential suppression pattern explains gene-level divergence (Jaccard = 0.06) within the framework of pathway-level convergence (r = 0.84). Evidence from TF Activity analysis (T4_S1).

Cell-Type Proportions by Condition

Estimated mean proportions of 21 cell types across four conditions. Lesional skin (AL) at baseline shows higher immune cell fractions (monocytes, neutrophils, macrophages, DCs) relative to non-lesional (AN). After treatment (m3), proportions shift toward non-lesional patterns as myeloid cells decrease and structural cells relatively increase.

21 cell types. Diff = activity difference (AL − AN); * = FDR < 0.05. Treatment changes are paired (m3 − m0) in AL skin. Dupilumab n=21, Cyclosporine n=8.

Immune Landscape

Myeloid dominance in AD lesions: Monocytes, neutrophils, macrophages, and dendritic cells are the most significantly elevated cell types in lesional skin, confirming the T4_S1 finding that SPI1/PU.1 (myeloid master regulator) is among the most differentially active TFs. Eosinophils (+0.332) and NK cells (+0.334) show the largest absolute differences, consistent with type 2 immune activation. Both treatments reduce myeloid infiltration (8 significant for dupilumab) and relatively increase structural cell proportions (fibroblasts, endothelial cells, pericytes), with treatment-level correlation r = 0.714 fitting the convergence hierarchy: gene 0.56 < cell-type 0.71 < TF 0.77 < pathway 0.84. Dupilumab paradoxically increases basophil scores (+0.042, FDR = 0.020), possibly reflecting compensatory basophil accumulation without activation after IL-4Rα blockade.

Basophil activity scores increase after dupilumab treatment (+0.042, FDR = 0.020) despite dupilumab blocking IL-4Rα which mediates basophil activation. This counter-intuitive finding may reflect compensatory basophil accumulation without activation — or a methodological limitation of marker-based deconvolution in distinguishing basophils from other granulocytes. Flow cytometry of basophils in dupilumab-treated AD patients would distinguish between increased numbers versus altered activation state. Evidence from Cell-Type Deconvolution analysis (T4_S2).

Module–Trait Correlations

Pearson correlation between module eigengenes and sample traits. * p < 0.05, ** p < 0.01, *** p < 0.001. ME1 (r = +0.63) and ME4 (r = -0.68) are the dominant disease-associated modules. Red = positive correlation with trait, Blue = negative.

Module Architecture

Two modules dominate the disease axis: M1 (1,417 genes, r = +0.63 with AL) captures the inflammatory/proliferative program — containing 94% of IFN-gamma response genes and 100% of E2F targets — with TYMP (kME = 0.913) as hub gene, confirming the original paper. M4 (392 genes, r = -0.68) is the barrier/structural module with BTC (kME = 0.907) as hub. Only ME1 is significantly reduced by dupilumab (Δ = -14.0, FDR = 0.020). Module-level convergence (r = 0.945) is the highest across all organizational levels: gene 0.56 < cell-type 0.71 < TF 0.77 < pathway 0.84 < per-sample pathway 0.90 < module 0.95.

M1 contains 94% of IFN-gamma response genes (45/48) and 100% of E2F targets (25/25), demonstrating that inflammation and proliferation are tightly co-regulated in AD lesional skin rather than independent transcriptional programs. This implies that anti-inflammatory treatment inherently normalizes keratinocyte proliferation, and vice versa — consistent with the dupilumab finding that anti-proliferative pathway suppression (E2F NES = −2.51) matches immune suppression in magnitude. Evidence from WGCNA (T5_S1) and Module Enrichment (T5_S3).

TYMP (thymidine phosphorylase, kME = 0.913) independently emerges as the #1 hub gene of the 1,417-gene inflammatory module M1, confirming the original Möbus et al. finding. TYMP catalyzes thymidine to thymine and deoxyribose-1-phosphate, with its product promoting angiogenesis and inflammatory cell chemotaxis. Its position at the center of the co-expression network makes it a potential biomarker and therapeutic target for the AD inflammatory program. Evidence from WGCNA (T5_S1–S2).

Biomarker Prediction: Actual vs Predicted Treatment Response

LOO-CV predicted vs actual composite inflammatory pathway response for 29 patients (13-gene LASSO panel). Pearson r = 0.428, R² = 0.113[95% CI: -0.456, 0.428]. Dashed line = perfect prediction. Points colored by treatment arm.

Positive coefficients (red) indicate higher baseline expression predicts greater inflammatory reduction. Top features: CLEC4OP (+0.23), IL19 (−0.21), LINC02610 (+0.19). CD274 (PD-L1) with negative coefficient suggests immune checkpoint status may modulate treatment response.

Unsupervised clustering reveals that 34% of AD patients (C1, n=19) harbor a low-inflammation molecular subtype in lesional skin, with significantly lower TNFa signaling, NF-kB activity, and STAT3 scores despite clinical lesional status. This suggests that AD molecular severity forms a continuum rather than a single disease state. The low cross-modality agreement (max ARI = 0.277) implies that subtypes manifest differently at gene, pathway, and regulatory levels — patients who cluster together by pathway activity do not necessarily group by TF activity. Evidence from Patient Clustering (T6_S1). Silhouette score = 0.131 and bootstrap ARI = 0.499 indicate borderline stability, consistent with biological continuity rather than discrete subtypes.

Biomarker Caveat: Exploratory Analysis

The 13-gene biomarker panel should be considered exploratory and hypothesis-generating. With only 29 patients, the LOO-CV R² of 0.113 has a 95% CI of [-0.456, 0.428] — crossing zero. Correlation pre-screening was performed outside the LOO loop, introducing optimistic bias. The observed correlation (r = 0.428) may partly reflect regression to the mean: patients with higher baseline inflammation have more “room to improve.” Independent validation in external cohorts is required before any clinical interpretation. The original paper’s candidate genes (IL4R, IL22, IL36A, S100A8, S100A9, IL13) show significant negative correlations with composite response, confirming this severity-response relationship.

Safety Gene Panel Heatmap

Heatmap of log₂FC values for 112 safety-relevant genes across 4 treatment comparisons. Rows grouped by safety panel. Red indicates upregulation, blue indicates downregulation. Stars (★) mark significant changes (FDR < 0.05, |log₂FC| > 0.5). The fibrosis panel shows striking cyclosporine-specific collagen upregulation (dark red cells), particularly COL1A1, COL1A2, and COL3A1 with log₂FC of 1.9–3.0.

35 flagged gene-comparison pairs (FDR < 0.05 and |log₂FC| > 0.5). Cyclosporine accounts for 29 of 35 flags.

Important Caveat

Transcriptomic changes do NOT equate to clinical toxicity. Gene expression shifts indicate pathway activation/suppression but require functional validation for safety assessment. These panels are designed to flag genes for further investigation, not to make definitive safety claims. Gene expression changes represent pathway-level signals that require protein-level, functional, and clinical validation before informing safety decisions.

Cyclosporine induces massive upregulation of collagens (COL1A1/COL1A2/COL3A1 log₂FC 1.9–3.0) plus LOX, LOXL2, FN1, and TGFB2, constituting 10 of its 22 flagged genes. This fibrotic signature is stronger in non-lesional than lesional skin (AN mean log₂FC = +0.75 vs AL +0.20), suggesting a direct cyclosporine effect on fibroblasts rather than secondary inflammation resolution. While cyclosporine-induced renal fibrosis is well documented, this skin-specific fibrotic program has not been previously characterized in AD transcriptomic studies. Evidence from Safety Panel Analysis (T7_S1).

Safety Profile Comparison

Dupilumab’s clean profile (5 flags) is consistent with its targeted IL-4Rα mechanism. The primary signal — GZMB suppression (log₂FC = −1.37) — represents normalized cytotoxic T-cell activation rather than unsafe immunosuppression. Cyclosporine’s broader footprint (22 flags) reflects its non-specific calcineurin inhibition affecting multiple pathways beyond target T-cell suppression. Notably, cyclosporine’s DNA damage response panel shows near-zero mean log₂FC (−0.09) with only 1/18 genes significant — consistent with literature showing its skin cancer risk operates through post-transcriptional mechanisms (XPC protein degradation, CypA-mediated checkpoint disruption) rather than transcriptional suppression. Both treatments suppress DUOX1 in lesional skin, a therapeutically beneficial effect: DUOX1 amplifies Th2 inflammation via a ROS-STAT6 positive feedback loop.

Treatment Convergence Hierarchy

Dupilumab vs cyclosporine correlation at each organizational level. The monotonic increase from gene (r = 0.56) to module (r = 0.95) reveals that different gene targets converge into shared biological programs at higher organizational levels.

Click any row to expand its full description. 7 high confidence, 5 medium confidence findings across 8 analytical tracks.

Mechanistic Integration

AD skin pathobiology is organized hierarchically: individual gene changes are highly treatment-specific (6% overlap), but these gene-level changes aggregate into remarkably convergent effects at the pathway, TF, and module levels (up to r = 0.95). This has profound implications: (1) gene-level biomarkers will be treatment-specific while pathway-level biomarkers will be treatment-general; (2) the 15-gene minimal shared signature (S100A12, CCL20, IL36A, SPRR2 family) represents the irreducible therapeutic core; (3) cyclosporine achieves equivalent pathway normalization through coordinated sub-threshold changes rather than large individual gene effects (9.1x gene-to-pathway amplification vs 3.4x for dupilumab). Safety analysis reveals a striking asymmetry: cyclosporine triggers 22 safety gene flags vs 5 for dupilumab, with a comprehensive pro-fibrotic signature stronger in non-lesional than lesional skin.

Testable Prediction

If the convergence hierarchy reflects genuine biology, then in an independent AD cohort treated with a third agent (e.g., JAK inhibitor), we should observe similarly low gene-level overlap with dupilumab/cyclosporine but high pathway-level correlation (r > 0.7), and the 15-gene shared signature should be reversed. This prediction is testable with publicly available transcriptomic data from baricitinib or upadacitinib trials in AD.

Study Limitations

22 analytical steps across 8 tracks. All steps used linear mixed-effects models (LMM) with patient-level random intercepts where applicable. FDR correction via Benjamini–Hochberg across all comparisons.