Analysis Report: DEMO

⚙️ Pipeline Logic (src/mapping_pipeline.py):

1. Resource Setup: Downloads the reference transcriptome from GENCODE, pulls the Salmon Apptainer image, and builds a Salmon index (if not already present).
2. Sample Discovery: Scans fastq_data/ for paired-end _R1.fastq.gz and matching _R2.fastq.gz sequence files, grouping them as distinct biological replicates.
3. Execution Guard (tmux): Spawns a background tmux session (rnaseq_<organism>) to protect mapping against potential network disconnects or process interruptions.
4. Salmon Quantification: Executes salmon quant within the container, validating mappings and correcting for GC bias (--validateMappings --gcBias).
5. Aggregation: Collects individual transcript count and TPM values from Salmon quant.sf output directories into combined matrices.
6. QC Visualizations: Generates sample PCA plot (PC1 vs PC2) on log2-transformed TPM data, and detected gene count barplot (genes with TPM > 1).

⚙️ Pipeline Logic (src/utils.py & src/mapping_pipeline.py):

• ID Parsing (parse_gene_ids): Parses Gencode pipe-delimited header rows (e.g., ENST|ENSG|...|Symbol|...) to map raw indices to Ensembl Gene IDs, retrieve gene symbols, and capture biotypes.
• Isoform Aggregation (aggregate_counts, load_and_aggregate_counts): Sums transcript-level raw counts or TPMs grouping them by Ensembl Gene ID or mapped Gene Symbol to obtain gene-level expression tables.

⚙️ Pipeline Logic (src/run_dge.py):

1. Isoform Aggregation: Extracts Ensembl Gene IDs from pipe-delimited raw row headers and aggregates isoform transcript counts to gene level by summing.
2. Gene Pruning: Filters out low-abundance genes with a total count sum ≤ 10 across all samples of the contrast.
3. Two-Pass Annotation: Employs GeneAnnotator (combining local caching, fast MyGene bulk queries, and concurrent Ensembl REST fallback) to fetch gene symbols, biotypes, and functional descriptions.
4. Statistical Fitting: Initialises DeseqDataSet (design='~Condition') and runs the DESeq2 pipeline (normalizes counts, estimates dispersion, fits negative binomial GLM, and refits outliers using Cook's distance).
5. Wald Test: Executes DeseqStats to perform pairwise Wald tests, computing Log2 Fold Changes, raw p-values, and BH-adjusted p-values (padj).
6. Export & Volcano Plots: Saves results to DGE CSVs. Generates 1x2 volcano plots showing all genes and significant genes (padj ≤ 0.05). Falls back to Log2FC histograms if replicates are insufficient ($N < 3$).

⚙️ Pipeline Logic for Unified Sample Comparison (src/dose_response_clustering.py):

To analyze how genes adapt across treatment doses without isolating the cell lines, we utilize a unified 8-dimensional clustering space. The mathematical workflow is as follows:

1. Compilation of Union DEG Set:
   First, we identify significant differentially expressed genes (padj ≤ 0.05) in any contrast of Cell Line A or Cell Line B. Taking the union of these genes ensures that any transcript exhibiting response dynamics in at least one cellular background is included.
2. Independent Baseline Standardization (Row Z-scoring):
   To isolate the relative shape (trend) of the response curve from static baseline abundance differences (since one cell line may naturally express a gene at a much higher absolute level than the other), the 4-point dose profile (0mM, 2mM, 5mM, 10mM) is averaged across replicates and row z-scored independently for each cell line: $$\mu_{\text{line}} = \frac{1}{4}\sum_{d} TPM_{line, d},\quad \sigma_{\text{line}} = \text{std}(TPM_{line})$$ $$Z_{line, d} = \frac{TPM_{line, d} - \mu_{\text{line}}}{\sigma_{\text{line}}}$$ This maps both profiles to a comparable scale with a mean of 0 and standard deviation of 1.
3. 8-Dimensional Trajectory Concatenation:
   For each gene, we concatenate the 4-point standardized profile of Cell Line A and the 4-point standardized profile of Cell Line B into a single 8-dimensional feature vector: $$\mathbf{V}_{gene} = [Z_{A,0}, Z_{A,2}, Z_{A,5}, Z_{A,10}, Z_{B,0}, Z_{B,2}, Z_{B,5}, Z_{B,10}]$$
4. Joint K-Means Clustering:
   K-Means clustering (using $K=8$ clusters) is performed directly on this 8-dimensional space. By clustering in 8D, the algorithm groups genes that share:
   • Shared Concordant Trends: E.g., genes that show a dose-dependent increase in both cell lines.
   • Divergent cell-line-specific Trends: E.g., genes that increase in Cell Line A but decrease or remain flat in Cell Line B.
5. Pathway Enrichment (Enrichr):
   The genes assigned to each cluster are submitted to the Enrichr API (gseapy.enrichr) to detect biological pathways (GO and KEGG) that correspond to these joint response patterns.

⚙️ Pipeline Logic (src/hdac_signature_analysis.py):

1. Signature Filtering: Subsets DGE results to predefined gene lists: HDACs (HDAC1-11, SIRT1-7), HATs (EP300, CREBBP, KATs), and canonical downstream target genes of HDAC inhibitors (e.g., CDKN1A, CCND1, metallothioneins).
2. Data Compilation: Merges log2FC, p-value, and adjusted p-value metrics for these genes across all contrasts and cell lines into a unified project dataframe.
3. Heatmap & Bar Plot Creation:
   • Generates a combined heatmap of HDAC/HAT Log2FC values across all cell lines (grouped by cell line and ordered by dose) to trace how epigenetic writers and erasers adapt across treatment concentrations.
   • Generates a combined bar plot showing dose-dependent fold changes of canonical HDACi target genes across all contrasts.
4. Targeted Volcano: Draws individual cell line volcano plots for their highest dose contrast (10mM) highlighting epigenetic regulators and targets.

⚙️ Pipeline Logic (src/mitochondrial_analysis.py):

1. OXPHOS-vs-Glycolysis Index:
   • Filters the TPM expression matrix for 11 core glycolysis genes and 30+ OXPHOS genes (Complexes I-V).
   • Computes the mean expression of OXPHOS and glycolysis genes for each sample.
   • Calculates the bioenergetic index as: OXPHOS_Mean / Glycolysis_Mean.
   • Merges all indexes into a unified project dataset and plots them side-by-side using different color fills for cell lines. Runs Mann-Whitney U test comparing the cell lines at each concentration.
2. Dynamics and Antioxidant Response:
   • Reads contrast-specific DGE results to retrieve Log2FC values for mitochondrial fusion genes (MFN1, MFN2, OPA1), fission genes (DNM1L, FIS1, MFF), and Nrf2 antioxidant response genes (SOD1, SOD2, CAT, GPX1, PRDX1, NFE2L2, HMOX1).
   • Generates unified project-level bar charts comparing these values across all sorted contrasts to profile mitochondrial structural dynamics and oxidative stress responses.

⚙️ Pipeline Logic (src/pathway_enrichment.py & src/export_pathway_json.py):

• GSEA Prerank (src/pathway_enrichment.py): Ranks all annotated genes in a contrast by log2FoldChange in descending order. Runs gseapy.prerank against reference databases (e.g. GO_Biological_Process_2025) using 1000 random permutations to calculate enrichment scores (ES), normalized enrichment scores (NES), and false discovery rate (FDR).
• GSVA (src/pathway_enrichment.py): Computes sample-level pathway enrichment scores from raw expression matrices using GSEAPy's GSVA algorithm to assess pathway activity variation across individual samples.
• Export Pathway JSON (src/export_pathway_json.py): Loads GSEA/GSVA report sheets, maps samples to experimental conditions using regex patterns, compiles parameters, and exports compressed payloads (gsea_data.json.gz, gsva_data.json.gz) for frontend display.

⚙️ Pipeline Logic (src/tf_metabolic_analysis.py & src/export_tf_metabolic_json.py):

1. Transcription Factor Activity (DoRothEA & CollecTRI):
   • Sample-level: Runs Multivariate Linear Model (dc.mt.mlm) on the sample-by-gene matrix to infer TF activity scores based on target gene expression.
   • Contrast-level: Runs Univariate Linear Model (dc.mt.ulm) using DGE contrast Log2FC vectors to determine TF activation state changes.
2. Network Graph Building (DoRothEA & CollecTRI):
   • Identifies contrasts in output folders and parses target link sheets (targets_<Contrast>.csv) to map source nodes (regulators) and target nodes (genes) with edge Mode of Regulation (mor).
   • Merges Ensembl descriptions and DGE results (Log2FC, padj) into gene nodes, and activity z-scores/p-values into regulatory nodes, exporting node-link schemas to compressed JSON files (tf_network.json.gz and collectri_network.json.gz).

⚙️ Pipeline Logic (src/tf_metabolic_analysis.py & src/export_tf_metabolic_json.py):

1. Signaling Pathway Profiling (PROGENy):
   • Sample-level: Runs dc.mt.mlm using PROGENy pathway weight matrices to determine pathway activation.
   • Contrast-level: Runs dc.mt.ulm on DGE contrast Log2FC vectors to determine signaling pathway activity changes.
2. Hallmark Pathways Activity:
   • Sample-level: Evaluates sample-specific pathway enrichment utilizing the AUCell algorithm (dc.mt.aucell).
   • Contrast-level: Runs dc.mt.ulm on DGE contrast Log2FC vectors.
3. Network Graph Building (PROGENy & Hallmark):
   • Parses target link sheets to map source nodes (pathways) and target genes, merging activity z-scores, p-values, and target LFC/padj metrics.
   • Exports compressed network schemas (progeny_network.json.gz and hallmark_network.json.gz) for visualization.