Skip to content

Instantly share code, notes, and snippets.

@wojtyniakAQ

wojtyniakAQ/notebook-ed17bd09-0916-40a1-8c1a-33bea019f2a1-paper-20260319-231703.ipynb

Created March 19, 2026 23:17
Show Gist options

  • Select an option

    Learn more about clone URLs
  • Save wojtyniakAQ/c9362720a75dc51e60039deeba464d20 to your computer and use it in GitHub Desktop.

Select an option

Learn more about clone URLs
Save wojtyniakAQ/c9362720a75dc51e60039deeba464d20 to your computer and use it in GitHub Desktop.
Download ZIP
Display the source blob
Display the rendered blob
Raw
notebook-ed17bd09-0916-40a1-8c1a-33bea019f2a1-paper-20260319-231703.ipynb
{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Developmental Regulation of Cell Type-Specific Transcription\n",
"\n",
"## Educational Notebook: RNA-seq Transcriptome Analysis Workflow\n",
"\n",
"**Paper**: *Developmental regulation of cell type-specific transcription by novel promoter-proximal sequence elements* \n",
"**Authors**: Dan Lu, Ho-Su Sin, Chenggang Lu, and Margaret T. Fuller \n",
"**Journal**: Genes & Development 34:663–677 (2020)\n",
"\n",
"---\n",
"\n",
"## Overview\n",
"\n",
"This notebook demonstrates the RNA-seq analysis workflow used in the paper to study transcriptional changes during the transition from spermatogonia to spermatocytes in *Drosophila*. The paper identified >3000 genes that are either newly expressed or expressed from alternative promoters in spermatocytes.\n",
"\n",
"### Key Findings:\n",
"- **1,155 down-regulated genes**: expressed ≥2-fold less in spermatocytes\n",
"- **1,841 off-to-on genes**: expressed >16-fold higher in spermatocytes\n",
"- **1,230 alternative promoter genes**: use different transcription start sites (TSS)\n",
"\n",
"### Experimental Design:\n",
"- **bam−/−**: Mutant flies where spermatogonia continue to proliferate (control)\n",
"- **48 hrPHS**: 48 hours post-heat shock (early spermatocytes)\n",
"- **72 hrPHS**: 72 hours post-heat shock (mature spermatocytes)\n",
"\n",
"---\n",
"\n",
"## Resource Constraints Note\n",
"\n",
"⚠️ **This is an educational demonstration** using small-scale synthetic data that runs in ~5-10 minutes. \n",
"For full-scale analysis with real data, you would need:\n",
"- High-performance computing cluster\n",
"- 16+ GB RAM\n",
"- Several hours of compute time\n",
"- GPU for some analyses (optional)\n",
"\n",
"**Note**: The paper used R/DESeq2, but this notebook uses Python/PyDESeq2 due to environment constraints. The methodology and results are equivalent."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 1. Install and Load Required Packages\n",
"\n",
"First, we'll install all necessary Python packages for the analysis."
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"\u001b[37m⠋\u001b[0m \u001b[2mResolving dependencies... \u001b[0m\r",
"\u001b[2K\u001b[37m⠙\u001b[0m \u001b[2mResolving dependencies... \u001b[0m\r",
"\u001b[2K\u001b[37m⠋\u001b[0m \u001b[2mResolving dependencies... \u001b[0m\r",
"\u001b[2K\u001b[37m⠙\u001b[0m \u001b[2mResolving dependencies... \u001b[0m"
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"\r",
"\u001b[2K\u001b[37m⠙\u001b[0m \u001b[2mpydeseq2==0.5.4 \u001b[0m"
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"\r",
"\u001b[2K\u001b[37m⠙\u001b[0m \u001b[2mnumpy==2.4.2 \u001b[0m\r",
"\u001b[2K\u001b[37m⠙\u001b[0m \u001b[2mpandas==3.0.1 \u001b[0m\r",
"\u001b[2K\u001b[37m⠙\u001b[0m \u001b[2mnumpy==2.4.2 \u001b[0m\r",
"\u001b[2K\u001b[37m⠙\u001b[0m \u001b[2mmatplotlib==3.10.8 \u001b[0m\r",
"\u001b[2K\u001b[37m⠙\u001b[0m \u001b[2mseaborn==0.13.2 \u001b[0m\r",
"\u001b[2K\u001b[37m⠙\u001b[0m \u001b[2mscipy==1.17.1 \u001b[0m\r",
"\u001b[2K\u001b[37m⠙\u001b[0m \u001b[2mscikit-learn==1.8.0 \u001b[0m\r",
"\u001b[2K\u001b[37m⠙\u001b[0m \u001b[2manndata==0.12.10 \u001b[0m\r",
"\u001b[2K\u001b[37m⠙\u001b[0m \u001b[2mnumpy==2.4.2 \u001b[0m\r",
"\u001b[2K\u001b[37m⠙\u001b[0m \u001b[2mmatplotlib==3.10.8 \u001b[0m\r",
"\u001b[2K\u001b[37m⠙\u001b[0m \u001b[2mseaborn==0.13.2 \u001b[0m\r",
"\u001b[2K\u001b[37m⠙\u001b[0m \u001b[2mscipy==1.17.1 \u001b[0m\r",
"\u001b[2K\u001b[37m⠙\u001b[0m \u001b[2mscikit-learn==1.8.0 \u001b[0m\r",
"\u001b[2K\u001b[37m⠙\u001b[0m \u001b[2manndata==0.12.9 \u001b[0m\r",
"\u001b[2K\u001b[37m⠙\u001b[0m \u001b[2manndata==0.12.8 \u001b[0m"
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"\r",
"\u001b[2K\u001b[37m⠙\u001b[0m \u001b[2manndata==0.12.7 \u001b[0m\r",
"\u001b[2K\u001b[37m⠙\u001b[0m \u001b[2manndata==0.12.6 \u001b[0m\r",
"\u001b[2K\u001b[37m⠙\u001b[0m \u001b[2mformulaic-contrasts==1.0.0 \u001b[0m\r",
"\u001b[2K\u001b[37m⠙\u001b[0m \u001b[2mformulaic==1.2.1 \u001b[0m\r",
"\u001b[2K\u001b[37m⠙\u001b[0m \u001b[2mpython-dateutil==2.9.0.post0 \u001b[0m"
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"\r",
"\u001b[2K\u001b[37m⠹\u001b[0m \u001b[2marray-api-compat==1.14.0 \u001b[0m"
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"\r",
"\u001b[2K\u001b[37m⠹\u001b[0m \u001b[2mwrapt==2.1.2 \u001b[0m\r",
"\u001b[2K\u001b[2mResolved \u001b[1m36 packages\u001b[0m \u001b[2min 282ms\u001b[0m\u001b[0m\r\n",
"\u001b[37m⠋\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/0) \r",
"\u001b[2K\u001b[37m⠋\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16) \r",
"\u001b[2K\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16) "
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"\r",
"\u001b[2K\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2mdonfig \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/21.09 KiB \u001b[1A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2mdonfig \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/21.09 KiB \u001b[1A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2mlegacy-api-wrap \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/9.94 KiB\r\n",
"\u001b[2mdonfig \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/21.09 KiB \u001b[2A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[2A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2mlegacy-api-wrap \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/9.94 KiB\r\n",
"\u001b[2mdonfig \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/21.09 KiB \u001b[2A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[2A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2mlegacy-api-wrap \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/9.94 KiB\r\n",
"\u001b[2mdonfig \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/21.09 KiB \u001b[2A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[2A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2msession-info \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/8.91 KiB\r\n",
"\u001b[2mlegacy-api-wrap \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/9.94 KiB\r\n",
"\u001b[2mdonfig \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/21.09 KiB \u001b[3A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[3A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2msession-info \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/8.91 KiB\r\n",
"\u001b[2mlegacy-api-wrap \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/9.94 KiB\r\n",
"\u001b[2mdonfig \u001b[0m \u001b[32m---------------------\u001b[30m\u001b[2m---------\u001b[0m\u001b[0m 14.91 KiB/21.09 KiB \u001b[3A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[3A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2msession-info \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/8.91 KiB\r\n",
"\u001b[2mlegacy-api-wrap \u001b[0m \u001b[32m------------------------------\u001b[30m\u001b[2m\u001b[0m\u001b[0m 9.94 KiB/9.94 KiB\r\n",
"\u001b[2mdonfig \u001b[0m \u001b[32m---------------------\u001b[30m\u001b[2m---------\u001b[0m\u001b[0m 14.91 KiB/21.09 KiB \u001b[3A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[3A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2msession-info \u001b[0m \u001b[32m------------------------------\u001b[30m\u001b[2m\u001b[0m\u001b[0m 8.91 KiB/8.91 KiB\r\n",
"\u001b[2mlegacy-api-wrap \u001b[0m \u001b[32m------------------------------\u001b[30m\u001b[2m\u001b[0m\u001b[0m 9.94 KiB/9.94 KiB\r\n",
"\u001b[2mdonfig \u001b[0m \u001b[32m---------------------\u001b[30m\u001b[2m---------\u001b[0m\u001b[0m 14.91 KiB/21.09 KiB \u001b[3A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[3A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2msession-info \u001b[0m \u001b[32m------------------------------\u001b[30m\u001b[2m\u001b[0m\u001b[0m 8.91 KiB/8.91 KiB\r\n",
"\u001b[2mlegacy-api-wrap \u001b[0m \u001b[32m------------------------------\u001b[30m\u001b[2m\u001b[0m\u001b[0m 9.94 KiB/9.94 KiB\r\n",
"\u001b[2mdonfig \u001b[0m \u001b[32m---------------------\u001b[30m\u001b[2m---------\u001b[0m\u001b[0m 14.91 KiB/21.09 KiB \u001b[3A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[3A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2msession-info \u001b[0m \u001b[32m------------------------------\u001b[30m\u001b[2m\u001b[0m\u001b[0m 8.91 KiB/8.91 KiB\r\n",
"\u001b[2mlegacy-api-wrap \u001b[0m \u001b[32m------------------------------\u001b[30m\u001b[2m\u001b[0m\u001b[0m 9.94 KiB/9.94 KiB\r\n",
"\u001b[2mdonfig \u001b[0m \u001b[32m---------------------\u001b[30m\u001b[2m---------\u001b[0m\u001b[0m 14.91 KiB/21.09 KiB \u001b[3A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[3A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2mlegacy-api-wrap \u001b[0m \u001b[32m------------------------------\u001b[30m\u001b[2m\u001b[0m\u001b[0m 9.94 KiB/9.94 KiB\r\n",
"\u001b[2mdonfig \u001b[0m \u001b[32m---------------------\u001b[30m\u001b[2m---------\u001b[0m\u001b[0m 14.91 KiB/21.09 KiB \u001b[2A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[2A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2mdonfig \u001b[0m \u001b[32m---------------------\u001b[30m\u001b[2m---------\u001b[0m\u001b[0m 14.91 KiB/21.09 KiB \u001b[1A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2mdonfig \u001b[0m \u001b[32m---------------------\u001b[30m\u001b[2m---------\u001b[0m\u001b[0m 14.91 KiB/21.09 KiB \u001b[1A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2minterface-meta \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/14.51 KiB\r\n",
"\u001b[2mdonfig \u001b[0m \u001b[32m---------------------\u001b[30m\u001b[2m---------\u001b[0m\u001b[0m 14.91 KiB/21.09 KiB \u001b[2A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[2A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2minterface-meta \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/14.51 KiB\r\n",
"\u001b[2mdonfig \u001b[0m \u001b[32m---------------------\u001b[30m\u001b[2m---------\u001b[0m\u001b[0m 14.91 KiB/21.09 KiB \u001b[2A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[2A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2minterface-meta \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/14.51 KiB\r\n",
"\u001b[2mdonfig \u001b[0m \u001b[32m---------------------\u001b[30m\u001b[2m---------\u001b[0m\u001b[0m 14.91 KiB/21.09 KiB \u001b[2A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[2A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2minterface-meta \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/14.51 KiB\r\n",
"\u001b[2mdonfig \u001b[0m \u001b[32m---------------------\u001b[30m\u001b[2m---------\u001b[0m\u001b[0m 14.91 KiB/21.09 KiB\r\n",
"\u001b[2marray-api-compat \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/58.71 KiB \u001b[3A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[3A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2minterface-meta \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/14.51 KiB\r\n",
"\u001b[2mdonfig \u001b[0m \u001b[32m---------------------\u001b[30m\u001b[2m---------\u001b[0m\u001b[0m 14.91 KiB/21.09 KiB\r\n",
"\u001b[2marray-api-compat \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/58.71 KiB \u001b[3A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[3A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2mformulaic-contrasts \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/9.82 KiB\r\n",
"\u001b[2mdonfig \u001b[0m \u001b[32m---------------------\u001b[30m\u001b[2m---------\u001b[0m\u001b[0m 14.91 KiB/21.09 KiB\r\n",
"\u001b[2mnatsort \u001b[0m \u001b[32m-----------\u001b[30m\u001b[2m-------------------\u001b[0m\u001b[0m 14.93 KiB/37.37 KiB\r\n",
"\u001b[2marray-api-compat \u001b[0m \u001b[32m-------\u001b[30m\u001b[2m-----------------------\u001b[0m\u001b[0m 14.87 KiB/58.71 KiB\r\n",
"\u001b[2mnarwhals \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/434.44 KiB \u001b[5A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[5A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2mformulaic-contrasts \u001b[0m \u001b[32m------------------------------\u001b[30m\u001b[2m\u001b[0m\u001b[0m 9.82 KiB/9.82 KiB\r\n",
"\u001b[2mnatsort \u001b[0m \u001b[32m-----------\u001b[30m\u001b[2m-------------------\u001b[0m\u001b[0m 14.93 KiB/37.37 KiB\r\n",
"\u001b[2marray-api-compat \u001b[0m \u001b[32m-------\u001b[30m\u001b[2m-----------------------\u001b[0m\u001b[0m 14.87 KiB/58.71 KiB\r\n",
"\u001b[2mnarwhals \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/434.44 KiB \u001b[4A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[4A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2mnatsort \u001b[0m \u001b[32m-----------\u001b[30m\u001b[2m-------------------\u001b[0m\u001b[0m 14.93 KiB/37.37 KiB\r\n",
"\u001b[2marray-api-compat \u001b[0m \u001b[32m-------\u001b[30m\u001b[2m-----------------------\u001b[0m\u001b[0m 14.87 KiB/58.71 KiB\r\n",
"\u001b[2mnarwhals \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/434.44 KiB \u001b[3A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[3A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2mnatsort \u001b[0m \u001b[32m------------------------------\u001b[30m\u001b[2m\u001b[0m\u001b[0m 37.37 KiB/37.37 KiB\r\n",
"\u001b[2mpydeseq2 \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/44.55 KiB\r\n",
"\u001b[2marray-api-compat \u001b[0m \u001b[32m---------------\u001b[30m\u001b[2m---------------\u001b[0m\u001b[0m 30.87 KiB/58.71 KiB\r\n",
"\u001b[2mstdlib-list \u001b[0m \u001b[32m-----\u001b[30m\u001b[2m-------------------------\u001b[0m\u001b[0m 14.91 KiB/85.56 KiB\r\n",
"\u001b[2mformulaic \u001b[0m \u001b[32m---\u001b[30m\u001b[2m---------------------------\u001b[0m\u001b[0m 14.85 KiB/114.54 KiB\r\n",
"\u001b[2mwrapt \u001b[0m \u001b[32m---\u001b[30m\u001b[2m---------------------------\u001b[0m\u001b[0m 14.92 KiB/118.60 KiB\r\n",
"\u001b[2manndata \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/168.22 KiB\r\n",
"\u001b[2mzarr \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/277.41 KiB\r\n",
"\u001b[2mnarwhals \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/434.44 KiB\r\n",
"\u001b[2mh5py \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 14.91 KiB/5.15 MiB\r\n",
"\u001b[2mnumcodecs \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/8.70 MiB \u001b[11A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[11A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2mpydeseq2 \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/44.55 KiB\r\n",
"\u001b[2marray-api-compat \u001b[0m \u001b[32m---------------\u001b[30m\u001b[2m---------------\u001b[0m\u001b[0m 30.87 KiB/58.71 KiB\r\n",
"\u001b[2mstdlib-list \u001b[0m \u001b[32m-----\u001b[30m\u001b[2m-------------------------\u001b[0m\u001b[0m 14.91 KiB/85.56 KiB\r\n",
"\u001b[2mformulaic \u001b[0m \u001b[32m---\u001b[30m\u001b[2m---------------------------\u001b[0m\u001b[0m 14.85 KiB/114.54 KiB\r\n",
"\u001b[2mwrapt \u001b[0m \u001b[32m---\u001b[30m\u001b[2m---------------------------\u001b[0m\u001b[0m 14.92 KiB/118.60 KiB\r\n",
"\u001b[2manndata \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/168.22 KiB\r\n",
"\u001b[2mzarr \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/277.41 KiB\r\n",
"\u001b[2mnarwhals \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/434.44 KiB\r\n",
"\u001b[2mh5py \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 14.91 KiB/5.15 MiB\r\n",
"\u001b[2mnumcodecs \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/8.70 MiB "
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"\u001b[10A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[10A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2mpydeseq2 \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/44.55 KiB\r\n",
"\u001b[2mstdlib-list \u001b[0m \u001b[32m-----\u001b[30m\u001b[2m-------------------------\u001b[0m\u001b[0m 14.91 KiB/85.56 KiB\r\n",
"\u001b[2mformulaic \u001b[0m \u001b[32m---\u001b[30m\u001b[2m---------------------------\u001b[0m\u001b[0m 14.85 KiB/114.54 KiB\r\n",
"\u001b[2mwrapt \u001b[0m \u001b[32m---\u001b[30m\u001b[2m---------------------------\u001b[0m\u001b[0m 14.92 KiB/118.60 KiB\r\n",
"\u001b[2manndata \u001b[0m \u001b[32m--\u001b[30m\u001b[2m----------------------------\u001b[0m\u001b[0m 14.91 KiB/168.22 KiB\r\n",
"\u001b[2mzarr \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 0 B/277.41 KiB\r\n",
"\u001b[2mnarwhals \u001b[0m \u001b[32m-\u001b[30m\u001b[2m-----------------------------\u001b[0m\u001b[0m 16.00 KiB/434.44 KiB\r\n",
"\u001b[2mh5py \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 14.91 KiB/5.15 MiB\r\n",
"\u001b[2mnumcodecs \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 14.91 KiB/8.70 MiB \u001b[9A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[9A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2mstdlib-list \u001b[0m \u001b[32m---------------------------\u001b[30m\u001b[2m---\u001b[0m\u001b[0m 78.91 KiB/85.56 KiB\r\n",
"\u001b[2mformulaic \u001b[0m \u001b[32m----------------\u001b[30m\u001b[2m--------------\u001b[0m\u001b[0m 62.85 KiB/114.54 KiB\r\n",
"\u001b[2mwrapt \u001b[0m \u001b[32m---------------\u001b[30m\u001b[2m---------------\u001b[0m\u001b[0m 62.92 KiB/118.60 KiB\r\n",
"\u001b[2manndata \u001b[0m \u001b[32m-----------\u001b[30m\u001b[2m-------------------\u001b[0m\u001b[0m 62.91 KiB/168.22 KiB\r\n",
"\u001b[2mzarr \u001b[0m \u001b[32m-----\u001b[30m\u001b[2m-------------------------\u001b[0m\u001b[0m 46.81 KiB/277.41 KiB\r\n",
"\u001b[2mnarwhals \u001b[0m \u001b[32m----\u001b[30m\u001b[2m--------------------------\u001b[0m\u001b[0m 62.67 KiB/434.44 KiB\r\n",
"\u001b[2mh5py \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 62.91 KiB/5.15 MiB\r\n",
"\u001b[2mnumcodecs \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 50.39 KiB/8.70 MiB \u001b[8A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[8A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2mformulaic \u001b[0m \u001b[32m------------------------\u001b[30m\u001b[2m------\u001b[0m\u001b[0m 94.85 KiB/114.54 KiB\r\n",
"\u001b[2mwrapt \u001b[0m \u001b[32m----------------------------\u001b[30m\u001b[2m--\u001b[0m\u001b[0m 110.92 KiB/118.60 KiB\r\n",
"\u001b[2manndata \u001b[0m \u001b[32m-------------------\u001b[30m\u001b[2m-----------\u001b[0m\u001b[0m 110.91 KiB/168.22 KiB\r\n",
"\u001b[2mzarr \u001b[0m \u001b[32m--------\u001b[30m\u001b[2m----------------------\u001b[0m\u001b[0m 78.92 KiB/277.41 KiB\r\n",
"\u001b[2mnarwhals \u001b[0m \u001b[32m-------\u001b[30m\u001b[2m-----------------------\u001b[0m\u001b[0m 110.67 KiB/434.44 KiB\r\n",
"\u001b[2mh5py \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 110.91 KiB/5.15 MiB\r\n",
"\u001b[2mnumcodecs \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 78.91 KiB/8.70 MiB \u001b[7A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[7A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2mformulaic \u001b[0m \u001b[32m-----------------------------\u001b[30m\u001b[2m-\u001b[0m\u001b[0m 111.83 KiB/114.54 KiB\r\n",
"\u001b[2manndata \u001b[0m \u001b[32m-------------------------\u001b[30m\u001b[2m-----\u001b[0m\u001b[0m 142.91 KiB/168.22 KiB\r\n",
"\u001b[2mzarr \u001b[0m \u001b[32m-----------\u001b[30m\u001b[2m-------------------\u001b[0m\u001b[0m 110.92 KiB/277.41 KiB\r\n",
"\u001b[2mnarwhals \u001b[0m \u001b[32m----------\u001b[30m\u001b[2m--------------------\u001b[0m\u001b[0m 158.67 KiB/434.44 KiB\r\n",
"\u001b[2mh5py \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 174.91 KiB/5.15 MiB\r\n",
"\u001b[2mnumcodecs \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 174.91 KiB/8.70 MiB \u001b[6A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[6A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2manndata \u001b[0m \u001b[32m-------------------------\u001b[30m\u001b[2m-----\u001b[0m\u001b[0m 142.91 KiB/168.22 KiB\r\n",
"\u001b[2mzarr \u001b[0m \u001b[32m-----------------\u001b[30m\u001b[2m-------------\u001b[0m\u001b[0m 158.92 KiB/277.41 KiB\r\n",
"\u001b[2mnarwhals \u001b[0m \u001b[32m----------\u001b[30m\u001b[2m--------------------\u001b[0m\u001b[0m 158.67 KiB/434.44 KiB\r\n",
"\u001b[2mh5py \u001b[0m \u001b[32m-\u001b[30m\u001b[2m-----------------------------\u001b[0m\u001b[0m 222.13 KiB/5.15 MiB\r\n",
"\u001b[2mnumcodecs \u001b[0m \u001b[32m\u001b[30m\u001b[2m------------------------------\u001b[0m\u001b[0m 220.57 KiB/8.70 MiB \u001b[5A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[5A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2manndata \u001b[0m \u001b[32m----------------------------\u001b[30m\u001b[2m--\u001b[0m\u001b[0m 158.91 KiB/168.22 KiB\r\n",
"\u001b[2mzarr \u001b[0m \u001b[32m------------------\u001b[30m\u001b[2m------------\u001b[0m\u001b[0m 174.92 KiB/277.41 KiB\r\n",
"\u001b[2mnarwhals \u001b[0m \u001b[32m-------------\u001b[30m\u001b[2m-----------------\u001b[0m\u001b[0m 190.67 KiB/434.44 KiB\r\n",
"\u001b[2mh5py \u001b[0m \u001b[32m---\u001b[30m\u001b[2m---------------------------\u001b[0m\u001b[0m 606.13 KiB/5.15 MiB\r\n",
"\u001b[2mnumcodecs \u001b[0m \u001b[32m--\u001b[30m\u001b[2m----------------------------\u001b[0m\u001b[0m 828.57 KiB/8.70 MiB "
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"\u001b[5A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[5A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2mzarr \u001b[0m \u001b[32m----------------------\u001b[30m\u001b[2m--------\u001b[0m\u001b[0m 206.92 KiB/277.41 KiB\r\n",
"\u001b[2mnarwhals \u001b[0m \u001b[32m--------------\u001b[30m\u001b[2m----------------\u001b[0m\u001b[0m 206.67 KiB/434.44 KiB\r\n",
"\u001b[2mh5py \u001b[0m \u001b[32m------\u001b[30m\u001b[2m------------------------\u001b[0m\u001b[0m 1.15 MiB/5.15 MiB\r\n",
"\u001b[2mnumcodecs \u001b[0m \u001b[32m-------------\u001b[30m\u001b[2m-----------------\u001b[0m\u001b[0m 3.80 MiB/8.70 MiB \u001b[4A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[4A\u001b[37m⠙\u001b[0m \u001b[2mPreparing packages...\u001b[0m (0/16)\r\n",
"\u001b[2mzarr \u001b[0m \u001b[32m----------------------\u001b[30m\u001b[2m--------\u001b[0m\u001b[0m 206.92 KiB/277.41 KiB\r\n",
"\u001b[2mnarwhals \u001b[0m \u001b[32m--------------\u001b[30m\u001b[2m----------------\u001b[0m\u001b[0m 206.67 KiB/434.44 KiB\r\n",
"\u001b[2mh5py \u001b[0m \u001b[32m--------\u001b[30m\u001b[2m----------------------\u001b[0m\u001b[0m 1.45 MiB/5.15 MiB\r\n",
"\u001b[2mnumcodecs \u001b[0m \u001b[32m-------------\u001b[30m\u001b[2m-----------------\u001b[0m\u001b[0m 3.87 MiB/8.70 MiB "
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"\u001b[4A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[4A\u001b[37m⠹\u001b[0m \u001b[2mPreparing packages...\u001b[0m (12/16)\r\n",
"\u001b[2mzarr \u001b[0m \u001b[32m-----------------------------\u001b[30m\u001b[2m-\u001b[0m\u001b[0m 270.92 KiB/277.41 KiB\r\n",
"\u001b[2mnarwhals \u001b[0m \u001b[32m------------------\u001b[30m\u001b[2m------------\u001b[0m\u001b[0m 270.67 KiB/434.44 KiB\r\n",
"\u001b[2mh5py \u001b[0m \u001b[32m---------------\u001b[30m\u001b[2m---------------\u001b[0m\u001b[0m 2.72 MiB/5.15 MiB\r\n",
"\u001b[2mnumcodecs \u001b[0m \u001b[32m---------------------\u001b[30m\u001b[2m---------\u001b[0m\u001b[0m 6.20 MiB/8.70 MiB "
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"\u001b[4A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[4A\u001b[37m⠹\u001b[0m \u001b[2mPreparing packages...\u001b[0m (12/16)\r\n",
"\u001b[2mnarwhals \u001b[0m \u001b[32m-------------------\u001b[30m\u001b[2m-----------\u001b[0m\u001b[0m 286.56 KiB/434.44 KiB\r\n",
"\u001b[2mh5py \u001b[0m \u001b[32m----------------\u001b[30m\u001b[2m--------------\u001b[0m\u001b[0m 2.75 MiB/5.15 MiB\r\n",
"\u001b[2mnumcodecs \u001b[0m \u001b[32m-----------------------------\u001b[30m\u001b[2m-\u001b[0m\u001b[0m 8.67 MiB/8.70 MiB \u001b[3A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[3A\u001b[37m⠹\u001b[0m \u001b[2mPreparing packages...\u001b[0m (12/16)\r\n",
"\u001b[2mnarwhals \u001b[0m \u001b[32m--------------------------\u001b[30m\u001b[2m----\u001b[0m\u001b[0m 382.67 KiB/434.44 KiB\r\n",
"\u001b[2mh5py \u001b[0m \u001b[32m----------------------\u001b[30m\u001b[2m--------\u001b[0m\u001b[0m 3.94 MiB/5.15 MiB\r\n",
"\u001b[2mnumcodecs \u001b[0m \u001b[32m-----------------------------\u001b[30m\u001b[2m-\u001b[0m\u001b[0m 8.70 MiB/8.70 MiB "
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"\u001b[3A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[3A\u001b[37m⠹\u001b[0m \u001b[2mPreparing packages...\u001b[0m (12/16)\r\n",
"\u001b[2mnarwhals \u001b[0m \u001b[32m---------------------------\u001b[30m\u001b[2m---\u001b[0m\u001b[0m 398.67 KiB/434.44 KiB\r\n",
"\u001b[2mh5py \u001b[0m \u001b[32m---------------------------\u001b[30m\u001b[2m---\u001b[0m\u001b[0m 4.78 MiB/5.15 MiB \u001b[2A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[2A\u001b[37m⠹\u001b[0m \u001b[2mPreparing packages...\u001b[0m (12/16)\r\n",
"\u001b[2mh5py \u001b[0m \u001b[32m---------------------------\u001b[30m\u001b[2m---\u001b[0m\u001b[0m 4.78 MiB/5.15 MiB \u001b[1A\r",
"\u001b[2K\u001b[1B\r",
"\u001b[2K\u001b[1A\u001b[37m⠹\u001b[0m \u001b[2mPreparing packages...\u001b[0m (12/16) \r",
"\u001b[2K\u001b[2mPrepared \u001b[1m16 packages\u001b[0m \u001b[2min 282ms\u001b[0m\u001b[0m\r\n",
"░░░░░░░░░░░░░░░░░░░░ [0/0] \u001b[2mInstalling wheels... \u001b[0m\r",
"\u001b[2K░░░░░░░░░░░░░░░░░░░░ [0/29] \u001b[2mInstalling wheels... \u001b[0m\r",
"\u001b[2K░░░░░░░░░░░░░░░░░░░░ [0/29] \u001b[2msession-info==1.0.1 \u001b[0m\r",
"\u001b[2K░░░░░░░░░░░░░░░░░░░░ [1/29] \u001b[2msession-info==1.0.1 \u001b[0m\r",
"\u001b[2K░░░░░░░░░░░░░░░░░░░░ [1/29] \u001b[2mpydeseq2==0.5.4 \u001b[0m\r",
"\u001b[2K█░░░░░░░░░░░░░░░░░░░ [2/29] \u001b[2mpydeseq2==0.5.4 \u001b[0m\r",
"\u001b[2K█░░░░░░░░░░░░░░░░░░░ [2/29] \u001b[2mlegacy-api-wrap==1.5 \u001b[0m\r",
"\u001b[2K██░░░░░░░░░░░░░░░░░░ [3/29] \u001b[2mlegacy-api-wrap==1.5 \u001b[0m\r",
"\u001b[2K██░░░░░░░░░░░░░░░░░░ [3/29] \u001b[2mstdlib-list==0.12.0 \u001b[0m\r",
"\u001b[2K██░░░░░░░░░░░░░░░░░░ [4/29] \u001b[2mstdlib-list==0.12.0 \u001b[0m\r",
"\u001b[2K██░░░░░░░░░░░░░░░░░░ [4/29] \u001b[2minterface-meta==1.3.0 \u001b[0m\r",
"\u001b[2K███░░░░░░░░░░░░░░░░░ [5/29] \u001b[2minterface-meta==1.3.0 \u001b[0m\r",
"\u001b[2K███░░░░░░░░░░░░░░░░░ [5/29] \u001b[2mwrapt==2.1.2 \u001b[0m\r",
"\u001b[2K████░░░░░░░░░░░░░░░░ [6/29] \u001b[2mwrapt==2.1.2 \u001b[0m\r",
"\u001b[2K████░░░░░░░░░░░░░░░░ [6/29] \u001b[2mdonfig==0.8.1.post1 \u001b[0m\r",
"\u001b[2K████░░░░░░░░░░░░░░░░ [7/29] \u001b[2mnarwhals==2.18.0 \u001b[0m\r",
"\u001b[2K█████░░░░░░░░░░░░░░░ [8/29] \u001b[2mnarwhals==2.18.0 \u001b[0m\r",
"\u001b[2K█████░░░░░░░░░░░░░░░ [8/29] \u001b[2mnarwhals==2.18.0 \u001b[0m\r",
"\u001b[2K█████░░░░░░░░░░░░░░░ [8/29] \u001b[2mformulaic-contrasts==1.0.0 \u001b[0m\r",
"\u001b[2K██████░░░░░░░░░░░░░░ [9/29] \u001b[2mformulaic-contrasts==1.0.0 \u001b[0m"
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"\r",
"\u001b[2K███████████████████░ [28/29] \u001b[2mpandas==3.0.1 \u001b[0m\r",
"\u001b[2K\u001b[2mInstalled \u001b[1m29 packages\u001b[0m \u001b[2min 67ms\u001b[0m\u001b[0m\r\n",
" \u001b[32m+\u001b[39m \u001b[1manndata\u001b[0m\u001b[2m==0.12.6\u001b[0m\r\n",
" \u001b[32m+\u001b[39m \u001b[1marray-api-compat\u001b[0m\u001b[2m==1.14.0\u001b[0m\r\n",
" \u001b[32m+\u001b[39m \u001b[1mcontourpy\u001b[0m\u001b[2m==1.3.3\u001b[0m\r\n",
" \u001b[32m+\u001b[39m \u001b[1mcycler\u001b[0m\u001b[2m==0.12.1\u001b[0m\r\n",
" \u001b[32m+\u001b[39m \u001b[1mdonfig\u001b[0m\u001b[2m==0.8.1.post1\u001b[0m\r\n",
" \u001b[32m+\u001b[39m \u001b[1mfonttools\u001b[0m\u001b[2m==4.62.1\u001b[0m\r\n",
" \u001b[32m+\u001b[39m \u001b[1mformulaic\u001b[0m\u001b[2m==1.2.1\u001b[0m\r\n",
" \u001b[32m+\u001b[39m \u001b[1mformulaic-contrasts\u001b[0m\u001b[2m==1.0.0\u001b[0m\r\n",
" \u001b[32m+\u001b[39m \u001b[1mh5py\u001b[0m\u001b[2m==3.16.0\u001b[0m\r\n",
" \u001b[32m+\u001b[39m \u001b[1minterface-meta\u001b[0m\u001b[2m==1.3.0\u001b[0m\r\n",
" \u001b[32m+\u001b[39m \u001b[1mjoblib\u001b[0m\u001b[2m==1.5.3\u001b[0m\r\n",
" \u001b[32m+\u001b[39m \u001b[1mkiwisolver\u001b[0m\u001b[2m==1.5.0\u001b[0m\r\n",
" \u001b[32m+\u001b[39m \u001b[1mlegacy-api-wrap\u001b[0m\u001b[2m==1.5\u001b[0m\r\n",
" \u001b[32m+\u001b[39m \u001b[1mmatplotlib\u001b[0m\u001b[2m==3.10.8\u001b[0m\r\n",
" \u001b[32m+\u001b[39m \u001b[1mnarwhals\u001b[0m\u001b[2m==2.18.0\u001b[0m\r\n",
" \u001b[32m+\u001b[39m \u001b[1mnatsort\u001b[0m\u001b[2m==8.4.0\u001b[0m\r\n",
" \u001b[32m+\u001b[39m \u001b[1mnumcodecs\u001b[0m\u001b[2m==0.16.5\u001b[0m\r\n",
" \u001b[32m+\u001b[39m \u001b[1mpandas\u001b[0m\u001b[2m==3.0.1\u001b[0m\r\n",
" \u001b[32m+\u001b[39m \u001b[1mpillow\u001b[0m\u001b[2m==12.1.1\u001b[0m\r\n",
" \u001b[32m+\u001b[39m \u001b[1mpydeseq2\u001b[0m\u001b[2m==0.5.4\u001b[0m\r\n",
" \u001b[32m+\u001b[39m \u001b[1mpyparsing\u001b[0m\u001b[2m==3.3.2\u001b[0m\r\n",
" \u001b[32m+\u001b[39m \u001b[1mscikit-learn\u001b[0m\u001b[2m==1.8.0\u001b[0m\r\n",
" \u001b[32m+\u001b[39m \u001b[1mscipy\u001b[0m\u001b[2m==1.17.1\u001b[0m\r\n",
" \u001b[32m+\u001b[39m \u001b[1mseaborn\u001b[0m\u001b[2m==0.13.2\u001b[0m\r\n",
" \u001b[32m+\u001b[39m \u001b[1msession-info\u001b[0m\u001b[2m==1.0.1\u001b[0m\r\n",
" \u001b[32m+\u001b[39m \u001b[1mstdlib-list\u001b[0m\u001b[2m==0.12.0\u001b[0m\r\n",
" \u001b[32m+\u001b[39m \u001b[1mthreadpoolctl\u001b[0m\u001b[2m==3.6.0\u001b[0m\r\n",
" \u001b[32m+\u001b[39m \u001b[1mwrapt\u001b[0m\u001b[2m==2.1.2\u001b[0m\r\n",
" \u001b[32m+\u001b[39m \u001b[1mzarr\u001b[0m\u001b[2m==3.1.5\u001b[0m\r\n"
]
}
],
"source": [
"# Install required packages\n",
"!uv pip install pydeseq2 numpy pandas matplotlib seaborn scipy scikit-learn anndata"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"✓ All packages loaded successfully!\n",
"NumPy version: 2.4.2\n",
"Pandas version: 3.0.1\n"
]
}
],
"source": [
"# Import all required libraries\n",
"import numpy as np\n",
"import pandas as pd\n",
"import matplotlib.pyplot as plt\n",
"import seaborn as sns\n",
"from scipy import stats\n",
"from sklearn.decomposition import PCA\n",
"from sklearn.preprocessing import StandardScaler\n",
"import warnings\n",
"\n",
"# Set random seed for reproducibility\n",
"np.random.seed(42)\n",
"\n",
"# Set plotting style\n",
"plt.style.use('seaborn-v0_8-darkgrid')\n",
"sns.set_palette(\"husl\")\n",
"warnings.filterwarnings('ignore')\n",
"\n",
"print(\"✓ All packages loaded successfully!\")\n",
"print(f\"NumPy version: {np.__version__}\")\n",
"print(f\"Pandas version: {pd.__version__}\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 2. Generate Synthetic RNA-seq Count Data\n",
"\n",
"We'll create synthetic RNA-seq count data that mimics the experimental design:\n",
"- **1,000 genes** (vs. ~9,371 in the real study)\n",
"- **2 replicates per condition** (vs. 2 in the real study)\n",
"- **3 conditions**: bam−/−, 48hrPHS, 72hrPHS\n",
"\n",
"The synthetic data will include:\n",
"- Down-regulated genes (~100)\n",
"- Off-to-on genes (~150)\n",
"- Unchanged genes (~750)\n",
"\n",
"### Scaling to Full Data:\n",
"For the full study, you would:\n",
"1. Download FASTQ files from GEO (accession GSE145975)\n",
"2. Align reads using STAR aligner (~2-4 hours)\n",
"3. Quantify with featureCounts or STAR (~30 min)\n",
"4. Process ~34-40 million reads per sample"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Generating synthetic RNA-seq count data...\n",
"Generated 1000 genes across 6 samples\n",
" - Down-regulated genes: 100\n",
" - Off-to-on genes: 150\n",
" - Unchanged genes: 750\n",
"\n",
"Count matrix shape: (1000, 6)\n",
"\n",
"Total counts per sample:\n",
"bam_rep1 880854\n",
"bam_rep2 881753\n",
"PHS48_rep1 906607\n",
"PHS48_rep2 919514\n",
"PHS72_rep1 986458\n",
"PHS72_rep2 965897\n",
"dtype: int64\n",
"\n",
"First few rows:\n",
" bam_rep1 bam_rep2 PHS48_rep1 PHS48_rep2 PHS72_rep1 PHS72_rep2\n",
"Gene_0000 958 1500 618 284 197 269\n",
"Gene_0001 2266 2487 2411 1587 543 511\n",
"Gene_0002 1014 1036 209 261 261 232\n",
"Gene_0003 766 668 314 284 181 198\n",
"Gene_0004 860 2717 763 688 335 367\n",
"\n",
"Sample information:\n",
" condition replicate\n",
"sample \n",
"bam_rep1 bam 1\n",
"bam_rep2 bam 2\n",
"PHS48_rep1 PHS48 1\n",
"PHS48_rep2 PHS48 2\n",
"PHS72_rep1 PHS72 1\n",
"PHS72_rep2 PHS72 2\n"
]
}
],
"source": [
"def generate_synthetic_counts(n_genes=1000, n_samples=6):\n",
" \"\"\"\n",
" Generate realistic RNA-seq count data mimicking the paper's design.\n",
" \n",
" Returns:\n",
" count_df: DataFrame with genes as rows and samples as columns\n",
" sample_info: DataFrame with sample metadata\n",
" \"\"\"\n",
" print(\"Generating synthetic RNA-seq count data...\")\n",
" \n",
" # Create gene names\n",
" gene_names = [f\"Gene_{i:04d}\" for i in range(n_genes)]\n",
" \n",
" # Define gene categories\n",
" n_down = 100 # Down-regulated genes\n",
" n_off_to_on = 150 # Off-to-on genes\n",
" n_unchanged = n_genes - n_down - n_off_to_on\n",
" \n",
" # Sample names\n",
" sample_names = ['bam_rep1', 'bam_rep2', 'PHS48_rep1', 'PHS48_rep2', 'PHS72_rep1', 'PHS72_rep2']\n",
" \n",
" # Initialize count matrix\n",
" counts = np.zeros((n_genes, n_samples), dtype=int)\n",
" \n",
" # Generate counts for each gene category\n",
" \n",
" # 1. Down-regulated genes (high in bam, low in spermatocytes)\n",
" for i in range(n_down):\n",
" base_mean = np.random.uniform(500, 3000)\n",
" # High in bam-/-\n",
" counts[i, 0:2] = np.random.negative_binomial(10, 10/(10+base_mean), 2)\n",
" # Low in 48hrPHS (2-4 fold down)\n",
" counts[i, 2:4] = np.random.negative_binomial(10, 10/(10+base_mean/np.random.uniform(2, 4)), 2)\n",
" # Very low in 72hrPHS (4-8 fold down)\n",
" counts[i, 4:6] = np.random.negative_binomial(10, 10/(10+base_mean/np.random.uniform(4, 8)), 2)\n",
" \n",
" # 2. Off-to-on genes (low in bam, high in spermatocytes)\n",
" for i in range(n_down, n_down + n_off_to_on):\n",
" base_mean = np.random.uniform(500, 3000)\n",
" # Very low in bam-/-\n",
" counts[i, 0:2] = np.random.negative_binomial(10, 10/(10+base_mean/20), 2)\n",
" # Moderate in 48hrPHS (8-16 fold up)\n",
" counts[i, 2:4] = np.random.negative_binomial(10, 10/(10+base_mean*np.random.uniform(0.4, 0.8)), 2)\n",
" # High in 72hrPHS (16-32 fold up from bam)\n",
" counts[i, 4:6] = np.random.negative_binomial(10, 10/(10+base_mean), 2)\n",
" \n",
" # 3. Unchanged genes (similar across conditions)\n",
" for i in range(n_down + n_off_to_on, n_genes):\n",
" base_mean = np.random.uniform(200, 2000)\n",
" counts[i, :] = np.random.negative_binomial(10, 10/(10+base_mean), n_samples)\n",
" \n",
" # Add some low-count genes to make it realistic\n",
" low_genes = np.random.choice(n_genes, 100, replace=False)\n",
" counts[low_genes, :] = np.random.negative_binomial(1, 1/(1+10), (100, n_samples))\n",
" \n",
" # Create count DataFrame\n",
" count_df = pd.DataFrame(counts, index=gene_names, columns=sample_names)\n",
" \n",
" # Create sample metadata\n",
" sample_info = pd.DataFrame({\n",
" 'sample': sample_names,\n",
" 'condition': ['bam', 'bam', 'PHS48', 'PHS48', 'PHS72', 'PHS72'],\n",
" 'replicate': [1, 2, 1, 2, 1, 2]\n",
" })\n",
" sample_info.set_index('sample', inplace=True)\n",
" \n",
" print(f\"Generated {n_genes} genes across {n_samples} samples\")\n",
" print(f\" - Down-regulated genes: {n_down}\")\n",
" print(f\" - Off-to-on genes: {n_off_to_on}\")\n",
" print(f\" - Unchanged genes: {n_unchanged}\")\n",
" \n",
" return count_df, sample_info\n",
"\n",
"# Generate synthetic data\n",
"count_matrix, sample_info = generate_synthetic_counts(n_genes=1000, n_samples=6)\n",
"\n",
"print(\"\\nCount matrix shape:\", count_matrix.shape)\n",
"print(\"\\nTotal counts per sample:\")\n",
"print(count_matrix.sum(axis=0))\n",
"print(\"\\nFirst few rows:\")\n",
"print(count_matrix.head())\n",
"print(\"\\nSample information:\")\n",
"print(sample_info)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 3. Quality Control: Explore Raw Count Data\n",
"\n",
"Before differential expression analysis, we'll examine the raw data quality."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# Calculate QC metrics\n",
"library_sizes = count_matrix.sum(axis=0)\n",
"genes_detected = (count_matrix > 0).sum(axis=0)\n",
"\n",
"qc_df = pd.DataFrame({\n",
" 'Sample': count_matrix.columns,\n",
" 'Condition': sample_info.loc[count_matrix.columns, 'condition'].values,\n",
" 'Library_Size': library_sizes.values,\n",
" 'Genes_Detected': genes_detected.values\n",
"})\n",
"\n",
"print(\"Quality Control Metrics:\")\n",
"print(qc_df)\n",
"\n",
"# Plot QC metrics\n",
"fig, axes = plt.subplots(1, 2, figsize=(14, 5))\n",
"\n",
"# Library sizes\n",
"colors = ['#66C2A5' if c == 'bam' else '#FC8D62' if c == 'PHS48' else '#8DA0CB' \n",
" for c in qc_df['Condition']]\n",
"axes[0].bar(qc_df['Sample'], qc_df['Library_Size']/1e6, color=colors)\n",
"axes[0].set_ylabel('Total Counts (millions)')\n",
"axes[0].set_title('Library Sizes (Total Counts)')\n",
"axes[0].tick_params(axis='x', rotation=45)\n",
"\n",
"# Genes detected\n",
"axes[1].bar(qc_df['Sample'], qc_df['Genes_Detected'], color=colors)\n",
"axes[1].set_ylabel('Genes with counts > 0')\n",
"axes[1].set_title('Number of Genes Detected')\n",
"axes[1].tick_params(axis='x', rotation=45)\n",
"\n",
"plt.tight_layout()\n",
"plt.show()\n",
"\n",
"print(\"\\n✓ QC plots generated successfully\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 4. Differential Expression Analysis with PyDESeq2\n",
"\n",
"Now we'll perform the main analysis using PyDESeq2 (Python implementation of DESeq2).\n",
"\n",
"### PyDESeq2 Workflow:\n",
"1. Filter low-count genes\n",
"2. Normalize counts and estimate size factors\n",
"3. Estimate dispersions\n",
"4. Perform statistical testing\n",
"5. Extract results for each comparison\n",
"\n",
"### Paper's Methods (Materials and Methods, page 674):\n",
"> \"Differential expression analyses were carried out using DESeq2\""
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"from pydeseq2.dds import DeseqDataSet\n",
"from pydeseq2.ds import DeseqStats\n",
"\n",
"print(\"Running PyDESeq2 differential expression analysis...\\n\")\n",
"\n",
"# Filter low-count genes (keep genes with at least 10 reads in at least 2 samples)\n",
"keep = (count_matrix >= 10).sum(axis=1) >= 2\n",
"count_filtered = count_matrix[keep]\n",
"\n",
"print(f\"Original genes: {len(count_matrix)}\")\n",
"print(f\"After filtering: {len(count_filtered)}\")\n",
"\n",
"# Prepare data for PyDESeq2\n",
"counts_df = count_filtered.copy()\n",
"\n",
"#Create metadata with proper format\n",
"metadata = sample_info.copy()\n",
"metadata['condition'] = pd.Categorical(metadata['condition'], \n",
" categories=['bam', 'PHS48', 'PHS72'])\n",
"\n",
"# Initialize DeseqDataSet\n",
"dds = DeseqDataSet(\n",
" counts=counts_df.T, # PyDESeq2 expects samples as rows\n",
" metadata=metadata,\n",
" design_factors='condition',\n",
" refit_cooks=True\n",
")\n",
"\n",
"print(\"\\nRunning DESeq2 analysis (this may take 1-2 minutes)...\")\n",
"dds.deseq2()\n",
"\n",
"print(\"\\n✓ DESeq2 analysis complete!\")\n",
"print(\"\\nSize factors (normalization):\")\n",
"print(dds.obsm['size_factors'])"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 5. Extract and Visualize Results\n",
"\n",
"Extract differential expression results for key comparisons."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# Extract results for 48hrPHS vs bam\n",
"print(\"Extracting differential expression results...\\n\")\n",
"\n",
"stat_res_48 = DeseqStats(dds, contrast=['condition', 'PHS48', 'bam'])\n",
"stat_res_48.summary()\n",
"results_48 = stat_res_48.results_df\n",
"\n",
"# Extract results for 72hrPHS vs bam\n",
"stat_res_72 = DeseqStats(dds, contrast=['condition', 'PHS72', 'bam'])\n",
"stat_res_72.summary()\n",
"results_72 = stat_res_72.results_df\n",
"\n",
"print(\"\\nResults for 48hrPHS vs bam:\")\n",
"print(f\"Total genes: {len(results_48)}\")\n",
"print(f\"Significant (padj < 0.05): {(results_48['padj'] < 0.05).sum()}\")\n",
"print(f\"Up-regulated: {((results_48['log2FoldChange'] > 0) & (results_48['padj'] < 0.05)).sum()}\")\n",
"print(f\"Down-regulated: {((results_48['log2FoldChange'] < 0) & (results_48['padj'] < 0.05)).sum()}\")\n",
"\n",
"print(\"\\nResults for 72hrPHS vs bam:\")\n",
"print(f\"Total genes: {len(results_72)}\")\n",
"print(f\"Significant (padj < 0.05): {(results_72['padj'] < 0.05).sum()}\")\n",
"print(f\"Up-regulated: {((results_72['log2FoldChange'] > 0) & (results_72['padj'] < 0.05)).sum()}\")\n",
"print(f\"Down-regulated: {((results_72['log2FoldChange'] < 0) & (results_72['padj'] < 0.05)).sum()}\")\n",
"\n",
"print(\"\\n✓ Results extracted successfully\")"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# Create MA plots (log ratio vs mean expression)\n",
"fig, axes = plt.subplots(1, 2, figsize=(14, 5))\n",
"\n",
"# MA plot for 48hrPHS vs bam\n",
"significant_48 = results_48['padj'] < 0.05\n",
"axes[0].scatter(results_48.loc[~significant_48, 'baseMean'], \n",
" results_48.loc[~significant_48, 'log2FoldChange'],\n",
" alpha=0.5, s=10, color='gray', label='Not significant')\n",
"axes[0].scatter(results_48.loc[significant_48, 'baseMean'], \n",
" results_48.loc[significant_48, 'log2FoldChange'],\n",
" alpha=0.7, s=10, color='blue', label='Significant')\n",
"axes[0].axhline(y=0, color='red', linestyle='--', linewidth=1)\n",
"axes[0].set_xlabel('Mean of normalized counts')\n",
"axes[0].set_ylabel('log2 Fold Change')\n",
"axes[0].set_title('MA Plot: 48hrPHS vs bam')\n",
"axes[0].set_xscale('log')\n",
"axes[0].set_ylim(-5, 5)\n",
"axes[0].legend()\n",
"axes[0].grid(True, alpha=0.3)\n",
"\n",
"# MA plot for 72hrPHS vs bam\n",
"significant_72 = results_72['padj'] < 0.05\n",
"axes[1].scatter(results_72.loc[~significant_72, 'baseMean'], \n",
" results_72.loc[~significant_72, 'log2FoldChange'],\n",
" alpha=0.5, s=10, color='gray', label='Not significant')\n",
"axes[1].scatter(results_72.loc[significant_72, 'baseMean'], \n",
" results_72.loc[significant_72, 'log2FoldChange'],\n",
" alpha=0.7, s=10, color='blue', label='Significant')\n",
"axes[1].axhline(y=0, color='red', linestyle='--', linewidth=1)\n",
"axes[1].set_xlabel('Mean of normalized counts')\n",
"axes[1].set_ylabel('log2 Fold Change')\n",
"axes[1].set_title('MA Plot: 72hrPHS vs bam')\n",
"axes[1].set_xscale('log')\n",
"axes[1].set_ylim(-5, 5)\n",
"axes[1].legend()\n",
"axes[1].grid(True, alpha=0.3)\n",
"\n",
"plt.tight_layout()\n",
"plt.show()\n",
"\n",
"print(\"\\n✓ MA plots generated\")\n",
"print(\" Blue points: significantly DE genes (padj < 0.05)\")\n",
"print(\" Gray points: non-significant genes\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 6. Gene Classification\n",
"\n",
"Following the paper's methodology, classify genes into categories:\n",
"\n",
"### From Results section (page 665):\n",
"1. **Down-regulated genes** (1,155 in paper): expressed ≥2-fold less in 48hrPHS or 72hrPHS vs bam−/−\n",
"2. **Off-to-on genes** (1,841 in paper): up-regulated >8-fold by 48hrPHS or >16-fold by 72hrPHS\n",
"3. **Alternative promoter genes** (1,230 in paper): expressed in both conditions but from different TSS\n",
"\n",
"Note: We can't classify alternative promoter genes without CAGE data, so we'll focus on categories 1 and 2."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"def classify_genes(results_48, results_72):\n",
" \"\"\"\n",
" Classify genes into expression categories based on fold change and significance.\n",
" \"\"\"\n",
" print(\"Classifying genes into expression categories...\\n\")\n",
" \n",
" # Get all genes\n",
" all_genes = results_72.index\n",
" \n",
" # Initialize classification\n",
" gene_class = pd.DataFrame({\n",
" 'gene': all_genes,\n",
" 'category': 'unchanged',\n",
" 'log2FC_48': results_48.loc[all_genes, 'log2FoldChange'].values,\n",
" 'log2FC_72': results_72.loc[all_genes, 'log2FoldChange'].values,\n",
" 'padj_48': results_48.loc[all_genes, 'padj'].values,\n",
" 'padj_72': results_72.loc[all_genes, 'padj'].values\n",
" }).set_index('gene')\n",
" \n",
" # Classify down-regulated genes (≥2-fold down, log2FC < -1)\n",
" down_48 = (gene_class['log2FC_48'] < -1) & (gene_class['padj_48'] < 0.05)\n",
" down_72 = (gene_class['log2FC_72'] < -1) & (gene_class['padj_72'] < 0.05)\n",
" down_genes = down_48 | down_72\n",
" gene_class.loc[down_genes, 'category'] = 'down-regulated'\n",
" \n",
" # Classify off-to-on genes\n",
" # >8-fold up by 48hrPHS (log2FC > 3) or >16-fold by 72hrPHS (log2FC > 4)\n",
" off_to_on_48 = (gene_class['log2FC_48'] > 3) & (gene_class['padj_48'] < 0.05)\n",
" off_to_on_72 = (gene_class['log2FC_72'] > 4) & (gene_class['padj_72'] < 0.05)\n",
" off_to_on_genes = off_to_on_48 | off_to_on_72\n",
" gene_class.loc[off_to_on_genes, 'category'] = 'off-to-on'\n",
" \n",
" # Summary\n",
" print(\"Gene classification summary:\")\n",
" print(f\" Total genes analyzed: {len(all_genes)}\")\n",
" print(f\" Down-regulated: {(gene_class['category'] == 'down-regulated').sum()}\")\n",
" print(f\" Off-to-on: {(gene_class['category'] == 'off-to-on').sum()}\")\n",
" print(f\" Unchanged: {(gene_class['category'] == 'unchanged').sum()}\")\n",
" \n",
" print(\"\\n✓ Gene classification complete\")\n",
" \n",
" return gene_class\n",
"\n",
"# Classify genes\n",
"gene_classes = classify_genes(results_48, results_72)\n",
"\n",
"# Display first few classified genes\n",
"print(\"\\nExample classified genes:\")\n",
"print(gene_classes.head(10))"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# Visualize gene categories\n",
"category_counts = gene_classes['category'].value_counts()\n",
"\n",
"fig, ax = plt.subplots(figsize=(8, 6))\n",
"colors_dict = {'down-regulated': '#E69F00', 'off-to-on': '#D55E00', 'unchanged': '#999999'}\n",
"colors = [colors_dict[cat] for cat in category_counts.index]\n",
"\n",
"bars = ax.bar(category_counts.index, category_counts.values, color=colors)\n",
"ax.set_ylabel('Number of Genes', fontsize=12)\n",
"ax.set_xlabel('Expression Category', fontsize=12)\n",
"ax.set_title('Gene Classification Summary', fontsize=14, fontweight='bold')\n",
"\n",
"# Add count labels on bars\n",
"for bar in bars:\n",
" height = bar.get_height()\n",
" ax.text(bar.get_x() + bar.get_width()/2., height,\n",
" f'{int(height)}',\n",
" ha='center', va='bottom', fontsize=11, fontweight='bold')\n",
"\n",
"plt.tight_layout()\n",
"plt.show()\n",
"\n",
"print(\"\\n✓ Category visualization complete\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 7. Scatter Plot Visualization\n",
"\n",
"Recreate Figure 1C from the paper showing RNA-seq expression levels in bam−/− vs 72hrPHS."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# Get normalized counts\n",
"normalized_counts = dds.layers['normed_counts'].T # Transpose back to genes x samples\n",
"normalized_df = pd.DataFrame(normalized_counts, \n",
" index=count_filtered.index, \n",
" columns=count_filtered.columns)\n",
"\n",
"# Calculate mean expression for each condition\n",
"mean_bam = normalized_df[['bam_rep1', 'bam_rep2']].mean(axis=1)\n",
"mean_72hrPHS = normalized_df[['PHS72_rep1', 'PHS72_rep2']].mean(axis=1)\n",
"\n",
"# Create scatter plot data\n",
"scatter_df = pd.DataFrame({\n",
" 'bam': np.log2(mean_bam + 1),\n",
" 'PHS72': np.log2(mean_72hrPHS + 1),\n",
" 'category': gene_classes.loc[mean_bam.index, 'category']\n",
"})\n",
"\n",
"# Create scatter plot\n",
"fig, ax = plt.subplots(figsize=(8, 8))\n",
"\n",
"colors_dict = {'down-regulated': '#E69F00', 'off-to-on': '#D55E00', 'unchanged': '#999999'}\n",
"\n",
"for category in ['unchanged', 'down-regulated', 'off-to-on']:\n",
" mask = scatter_df['category'] == category\n",
" ax.scatter(scatter_df.loc[mask, 'bam'], \n",
" scatter_df.loc[mask, 'PHS72'],\n",
" c=colors_dict[category], \n",
" alpha=0.5, \n",
" s=20, \n",
" label=category.replace('-', ' ').title())\n",
"\n",
"# Add diagonal line\n",
"max_val = max(scatter_df['bam'].max(), scatter_df['PHS72'].max())\n",
"ax.plot([0, max_val], [0, max_val], 'k--', linewidth=1, alpha=0.5)\n",
"\n",
"ax.set_xlabel('log2(normalized counts) in bam−/−', fontsize=12)\n",
"ax.set_ylabel('log2(normalized counts) in 72hrPHS', fontsize=12)\n",
"ax.set_title('RNA-seq Expression: bam−/− vs 72hrPHS\\n(Recreating Figure 1C)', \n",
" fontsize=14, fontweight='bold')\n",
"ax.legend(loc='upper left', frameon=True, fancybox=True)\n",
"ax.grid(True, alpha=0.3)\n",
"\n",
"plt.tight_layout()\n",
"plt.show()\n",
"\n",
"print(\"\\n✓ Scatter plot visualization complete\")\n",
"print(\"\\nInterpretation:\")\n",
"print(\" - Points on diagonal: similar expression in both conditions\")\n",
"print(\" - Points below diagonal (orange): down-regulated in spermatocytes\")\n",
"print(\" - Points above diagonal (red): off-to-on genes (newly expressed)\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 8. Heatmap of Top Differentially Expressed Genes\n",
"\n",
"Visualize expression patterns of the most significantly changed genes."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# Get top 50 most variable genes\n",
"gene_var = normalized_df.var(axis=1)\n",
"top_genes = gene_var.nlargest(50).index\n",
"\n",
"# Get expression matrix for top genes and scale\n",
"heatmap_data = normalized_df.loc[top_genes]\n",
"heatmap_scaled = (heatmap_data.T - heatmap_data.T.mean()) / heatmap_data.T.std()\n",
"\n",
"# Create heatmap\n",
"fig, ax = plt.subplots(figsize=(10, 12))\n",
"\n",
"# Create custom colormap\n",
"cmap = sns.diverging_palette(250, 10, as_cmap=True)\n",
"\n",
"# Plot heatmap\n",
"sns.heatmap(heatmap_scaled.T, \n",
" cmap=cmap, \n",
" center=0, \n",
" yticklabels=False,\n",
" cbar_kws={'label': 'Scaled expression'},\n",
" ax=ax)\n",
"\n",
"ax.set_title('Top 50 Variable Genes Across Conditions', fontsize=14, fontweight='bold')\n",
"ax.set_xlabel('Samples', fontsize=12)\n",
"ax.set_ylabel('Genes', fontsize=12)\n",
"\n",
"plt.tight_layout()\n",
"plt.show()\n",
"\n",
"print(\"\\n✓ Heatmap visualization complete\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 9. PCA Analysis\n",
"\n",
"Principal Component Analysis to visualize overall sample relationships."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# Perform PCA on normalized, log-transformed counts\n",
"log_counts = np.log2(normalized_df + 1)\n",
"\n",
"# Standardize\n",
"scaler = StandardScaler()\n",
"log_counts_scaled = scaler.fit_transform(log_counts.T) # Samples as rows\n",
"\n",
"# PCA\n",
"pca = PCA(n_components=2)\n",
"pca_result = pca.fit_transform(log_counts_scaled)\n",
"\n",
"# Create PCA DataFrame\n",
"pca_df = pd.DataFrame({\n",
" 'PC1': pca_result[:, 0],\n",
" 'PC2': pca_result[:, 1],\n",
" 'sample': count_filtered.columns,\n",
" 'condition': sample_info.loc[count_filtered.columns, 'condition'].values\n",
"})\n",
"\n",
"# Plot PCA\n",
"fig, ax = plt.subplots(figsize=(10, 8))\n",
"\n",
"colors_dict = {'bam': '#66C2A5', 'PHS48': '#FC8D62', 'PHS72': '#8DA0CB'}\n",
"labels_dict = {'bam': 'bam−/−', 'PHS48': '48hrPHS', 'PHS72': '72hrPHS'}\n",
"\n",
"for condition in ['bam', 'PHS48', 'PHS72']:\n",
" mask = pca_df['condition'] == condition\n",
" ax.scatter(pca_df.loc[mask, 'PC1'], \n",
" pca_df.loc[mask, 'PC2'],\n",
" c=colors_dict[condition], \n",
" s=200, \n",
" label=labels_dict[condition],\n",
" edgecolors='black',\n",
" linewidths=1.5)\n",
" \n",
" # Add sample labels\n",
" for idx, row in pca_df[mask].iterrows():\n",
" ax.annotate(row['sample'], \n",
" (row['PC1'], row['PC2']),\n",
" xytext=(5, 5), \n",
" textcoords='offset points',\n",
" fontsize=9)\n",
"\n",
"ax.set_xlabel(f'PC1 ({pca.explained_variance_ratio_[0]*100:.1f}% variance)', fontsize=12)\n",
"ax.set_ylabel(f'PC2 ({pca.explained_variance_ratio_[1]*100:.1f}% variance)', fontsize=12)\n",
"ax.set_title('PCA of RNA-seq Samples', fontsize=14, fontweight='bold')\n",
"ax.legend(fontsize=11, frameon=True, fancybox=True)\n",
"ax.grid(True, alpha=0.3)\n",
"\n",
"plt.tight_layout()\n",
"plt.show()\n",
"\n",
"print(\"\\n✓ PCA analysis complete\")\n",
"print(\"\\nInterpretation:\")\n",
"print(\" - PC1 separates conditions based on developmental stage\")\n",
"print(\" - Replicates cluster together, indicating good reproducibility\")\n",
"print(\" - Clear separation between spermatogonia and spermatocytes\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 10. Volcano Plots\n",
"\n",
"Visualize statistical significance vs fold change for each comparison."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# Prepare data for volcano plots\n",
"def create_volcano_data(results):\n",
" volcano_df = results.copy()\n",
" volcano_df['-log10(padj)'] = -np.log10(volcano_df['padj'] + 1e-300) # Add small value to avoid log(0)\n",
" \n",
" # Define significance categories\n",
" volcano_df['significance'] = 'Not significant'\n",
" volcano_df.loc[(volcano_df['padj'] < 0.05) & (volcano_df['log2FoldChange'] > 3), 'significance'] = 'Up (>8-fold)'\n",
" volcano_df.loc[(volcano_df['padj'] < 0.05) & (volcano_df['log2FoldChange'] < -1), 'significance'] = 'Down (>2-fold)'\n",
" volcano_df.loc[(volcano_df['padj'] < 0.05) & \n",
" (volcano_df['log2FoldChange'] >= -1) & \n",
" (volcano_df['log2FoldChange'] <= 3), 'significance'] = 'Significant'\n",
" \n",
" return volcano_df\n",
"\n",
"volcano_48 = create_volcano_data(results_48)\n",
"volcano_72 = create_volcano_data(results_72)\n",
"\n",
"# Create volcano plots\n",
"fig, axes = plt.subplots(1, 2, figsize=(16, 6))\n",
"\n",
"colors_dict = {\n",
" 'Up (>8-fold)': '#D55E00',\n",
" 'Down (>2-fold)': '#E69F00',\n",
" 'Significant': '#56B4E9',\n",
" 'Not significant': '#999999'\n",
"}\n",
"\n",
"# Plot 48hrPHS vs bam\n",
"for sig_type in ['Not significant', 'Significant', 'Down (>2-fold)', 'Up (>8-fold)']:\n",
" mask = volcano_48['significance'] == sig_type\n",
" axes[0].scatter(volcano_48.loc[mask, 'log2FoldChange'], \n",
" volcano_48.loc[mask, '-log10(padj)'],\n",
" c=colors_dict[sig_type], \n",
" alpha=0.6, \n",
" s=15, \n",
" label=sig_type)\n",
"\n",
"axes[0].axvline(x=-1, linestyle='--', color='gray', linewidth=1)\n",
"axes[0].axvline(x=3, linestyle='--', color='gray', linewidth=1)\n",
"axes[0].axhline(y=-np.log10(0.05), linestyle='--', color='gray', linewidth=1)\n",
"axes[0].set_xlabel('log2 Fold Change', fontsize=12)\n",
"axes[0].set_ylabel('-log10(adjusted p-value)', fontsize=12)\n",
"axes[0].set_title('48hrPHS vs bam−/−', fontsize=13, fontweight='bold')\n",
"axes[0].legend(fontsize=9, frameon=True)\n",
"axes[0].grid(True, alpha=0.3)\n",
"\n",
"# Plot 72hrPHS vs bam\n",
"for sig_type in ['Not significant', 'Significant', 'Down (>2-fold)', 'Up (>8-fold)']:\n",
" mask = volcano_72['significance'] == sig_type\n",
" axes[1].scatter(volcano_72.loc[mask, 'log2FoldChange'], \n",
" volcano_72.loc[mask, '-log10(padj)'],\n",
" c=colors_dict[sig_type], \n",
" alpha=0.6, \n",
" s=15, \n",
" label=sig_type)\n",
"\n",
"axes[1].axvline(x=-1, linestyle='--', color='gray', linewidth=1)\n",
"axes[1].axvline(x=3, linestyle='--', color='gray', linewidth=1)\n",
"axes[1].axhline(y=-np.log10(0.05), linestyle='--', color='gray', linewidth=1)\n",
"axes[1].set_xlabel('log2 Fold Change', fontsize=12)\n",
"axes[1].set_ylabel('-log10(adjusted p-value)', fontsize=12)\n",
"axes[1].set_title('72hrPHS vs bam−/−', fontsize=13, fontweight='bold')\n",
"axes[1].legend(fontsize=9, frameon=True)\n",
"axes[1].grid(True, alpha=0.3)\n",
"\n",
"plt.tight_layout()\n",
"plt.show()\n",
"\n",
"print(\"\\n✓ Volcano plots generated\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 11. Summary Statistics\n",
"\n",
"Generate comprehensive summary statistics."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# Create comprehensive results table\n",
"final_results = pd.DataFrame({\n",
" 'gene': gene_classes.index,\n",
" 'baseMean': results_72.loc[gene_classes.index, 'baseMean'].values,\n",
" 'log2FC_48vs_bam': gene_classes['log2FC_48'].values,\n",
" 'padj_48vs_bam': gene_classes['padj_48'].values,\n",
" 'log2FC_72vs_bam': gene_classes['log2FC_72'].values,\n",
" 'padj_72vs_bam': gene_classes['padj_72'].values,\n",
" 'category': gene_classes['category'].values\n",
"})\n",
"\n",
"print(\"=== FINAL SUMMARY STATISTICS ===\\n\")\n",
"\n",
"print(\"Gene Classification:\")\n",
"print(final_results['category'].value_counts())\n",
"\n",
"print(\"\\nTop 10 Down-regulated Genes (72hrPHS vs bam):\")\n",
"top_down = final_results[\n",
" final_results['category'] == 'down-regulated'\n",
"].sort_values('log2FC_72vs_bam').head(10)\n",
"print(top_down[['gene', 'log2FC_72vs_bam', 'padj_72vs_bam']].to_string(index=False))\n",
"\n",
"print(\"\\nTop 10 Off-to-on Genes (72hrPHS vs bam):\")\n",
"top_up = final_results[\n",
" final_results['category'] == 'off-to-on'\n",
"].sort_values('log2FC_72vs_bam', ascending=False).head(10)\n",
"print(top_up[['gene', 'log2FC_72vs_bam', 'padj_72vs_bam']].to_string(index=False))\n",
"\n",
"# Expression level summary\n",
"print(\"\\n\\nExpression Level Summary by Category:\")\n",
"expression_summary = final_results.groupby('category').agg({\n",
" 'gene': 'count',\n",
" 'baseMean': 'mean',\n",
" 'log2FC_72vs_bam': 'median'\n",
"}).rename(columns={'gene': 'n_genes', 'baseMean': 'mean_baseMean', \n",
" 'log2FC_72vs_bam': 'median_log2FC_72'})\n",
"print(expression_summary)\n",
"\n",
"print(\"\\n✓ Summary statistics complete\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 12. Analysis Complete!\n",
"\n",
"### Results Summary"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"print(\"=\"*60)\n",
"print(\"ANALYSIS COMPLETE\")\n",
"print(\"=\"*60)\n",
"print(\"\\nResults Summary:\")\n",
"print(f\" Total genes analyzed: {len(final_results)}\")\n",
"print(f\" Down-regulated: {(final_results['category'] == 'down-regulated').sum()}\")\n",
"print(f\" Off-to-on: {(final_results['category'] == 'off-to-on').sum()}\")\n",
"print(f\" Unchanged: {(final_results['category'] == 'unchanged').sum()}\")\n",
"\n",
"print(\"\\n✓ All analyses completed successfully!\")\n",
"print(\"\\nThis notebook demonstrated:\")\n",
"print(\" 1. Synthetic RNA-seq data generation\")\n",
"print(\" 2. Quality control analysis\")\n",
"print(\" 3. DESeq2 differential expression analysis\")\n",
"print(\" 4. Gene classification\")\n",
"print(\" 5. Multiple visualization techniques\")\n",
"print(\" 6. PCA for sample relationships\")\n",
"print(\"\\nThe workflow follows Lu et al. (2020) Genes & Development 34:663–677\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 13. Scaling to Full Experiments\n",
"\n",
"### How to Run This Analysis with Real Data\n",
"\n",
"This notebook demonstrated the workflow with small-scale synthetic data. To run with the actual paper data:\n",
"\n",
"#### 1. Data Acquisition\n",
"```bash\n",
"# Download data from GEO (requires ~50GB storage)\n",
"# Accession: GSE145975\n",
"wget ftp://ftp.ncbi.nlm.nih.gov/geo/series/GSE145nnn/GSE145975/suppl/GSE145975_RAW.tar\n",
"```\n",
"\n",
"#### 2. Read Alignment (requires ~4-8 hours, 16GB RAM)\n",
"```bash\n",
"# Using STAR aligner (as in paper)\n",
"STAR --runThreadN 8 \\\n",
" --genomeDir dm6_index \\\n",
" --readFilesIn sample_R1.fastq.gz sample_R2.fastq.gz \\\n",
" --readFilesCommand zcat \\\n",
" --outSAMtype BAM SortedByCoordinate \\\n",
" --quantMode GeneCounts\n",
"```\n",
"\n",
"#### 3. Count Quantification\n",
"```bash\n",
"# STAR outputs gene counts directly\n",
"# Or use featureCounts for more control\n",
"featureCounts -a genes.gtf -o counts.txt *.bam\n",
"```\n",
"\n",
"#### 4. DESeq2 Analysis (similar to this notebook)\n",
"```python\n",
"# Load real count data\n",
"count_data = pd.read_csv(\"counts.txt\", sep=\"\\t\", index_col=0)\n",
"# Continue with PyDESeq2 as shown above\n",
"```\n",
"\n",
"#### 5. Additional Analyses from the Paper\n",
"\n",
"The paper also performed:\n",
"- **CAGE-seq**: For precise TSS mapping (requires specialized library prep)\n",
"- **ATAC-seq**: For chromatin accessibility (requires ~10,000 cells per sample)\n",
"- **ChIP-seq**: For tMAC binding sites (requires antibodies)\n",
"- **Motif analysis**: Using MEME-ChIP suite\n",
"\n",
"### Computational Requirements for Full Analysis:\n",
"- **RAM**: 16-32 GB recommended\n",
"- **CPU**: 8+ cores\n",
"- **Storage**: 100+ GB\n",
"- **Time**: 1-2 days for complete pipeline\n",
"- **Software**: STAR, SAMtools, PyDESeq2, MEME suite"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 14. Key Biological Insights\n",
"\n",
"### Main Findings from the Paper:\n",
"\n",
"1. **Massive Transcriptional Reprogramming**: >3,000 genes (~30% of testis transcriptome) use new promoters during spermatocyte differentiation\n",
"\n",
"2. **Novel Promoter Elements**: Spermatocyte-specific promoters lack canonical TATA/DPE motifs but are enriched for:\n",
" - **tMAC-ChIP motif**: Binds the tMAC complex\n",
" - **Achi/Vis motif** (TGTCA): TALE-class homeodomain binding site\n",
" - **Inr motif** (TCA): At transcription start site\n",
" - **ACA motif**: At positions +26 to +30\n",
" - **CNAAATT motif**: At positions +29 to +60\n",
"\n",
"3. **tMAC Complex Function**: \n",
" - Opens chromatin at spermatocyte-specific promoters\n",
" - Creates ~100bp nucleosome-free region\n",
" - Required for most off-to-on gene expression\n",
"\n",
"4. **Promoter Architecture**: \n",
" - Upstream motifs (tMAC, Achi/Vis) recruit complexes and open chromatin\n",
" - Downstream motifs (ACA, CNAAATT) help specify efficient TSS usage\n",
" - Narrow-high promoters have optimal spacing between all motifs\n",
"\n",
"### Broader Significance:\n",
"- Demonstrates cell type-specific transcription can be driven by **promoter-proximal elements**, not just distal enhancers\n",
"- Shows how tissue-specific protein complexes (tMAC) can orchestrate developmental gene expression programs\n",
"- Reveals that alternative promoters are a major mechanism for stage-specific gene regulation"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## 15. Conclusion\n",
"\n",
"This notebook has demonstrated the core RNA-seq analysis workflow from Lu et al. (2020):\n",
"\n",
"✅ **Completed Steps:**\n",
"1. Generated synthetic RNA-seq count data mimicking the experimental design\n",
"2. Performed quality control on count data\n",
"3. Ran PyDESeq2 differential expression analysis\n",
"4. Classified genes into down-regulated, off-to-on, and unchanged categories\n",
"5. Visualized results with MA plots, scatter plots, volcano plots, and heatmaps\n",
"6. Performed PCA to assess sample relationships\n",
"\n",
"### Learning Outcomes:\n",
"- Understanding of developmental transcriptome analysis\n",
"- Experience with DESeq2-style differential expression analysis\n",
"- Gene classification strategies\n",
"- Data visualization for genomics\n",
"\n",
"### Next Steps for Full Analysis:\n",
"1. Obtain real RNA-seq data from GEO (GSE145975)\n",
"2. Perform read alignment with STAR\n",
"3. Run this workflow on real data\n",
"4. Add CAGE-seq analysis for TSS mapping\n",
"5. Add ATAC-seq for chromatin accessibility\n",
"6. Perform motif enrichment analysis\n",
"\n",
"---\n",
"\n",
"**References:**\n",
"- Lu et al. (2020) Genes & Development 34:663–677\n",
"- Love et al. (2014) Genome Biology 15:550 (DESeq2)\n",
"- GEO Accession: GSE145975"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.10.0"
}
},
"nbformat": 4,
"nbformat_minor": 4
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment