Quickstart + Usage¶

After installing Consenrich, check the build and try running from the command line with consenrich -h

Getting Started: H3K27ac ChIP-seq¶

This minimal example estimates a consensus from four ENCODE epidermis ChIP-seq experiments with matched input controls.

Input Data¶

H3K27ac demo inputs¶
Experiment	Biosample	H3K27ac alignment	Control alignment
ENCSR214UZE	Epidermis/Female/71	ENCFF793ZHL.bam	ENCFF444WVG.bam
ENCSR334DRN	Epidermis/Male/67	ENCFF647VPO.bam	ENCFF619NYP.bam
ENCSR340ZTB	Epidermis/Female/80	ENCFF809VKT.bam	ENCFF898LKJ.bam
ENCSR386CKJ	Epidermis/Male/75	ENCFF295EFL.bam	ENCFF490MWV.bam

Download Alignments¶

Copy and paste the following into your terminal to download and index the BAM files. You can use curl -L -O <URL> in place of wget <URL> if wget is not available.

encodeFiles=https://www.encodeproject.org/files
for file in ENCFF793ZHL ENCFF647VPO ENCFF809VKT ENCFF295EFL; do
    wget "$encodeFiles/$file/@@download/$file.bam"
done
for ctrl in ENCFF444WVG ENCFF619NYP ENCFF898LKJ ENCFF490MWV; do
    wget "$encodeFiles/$ctrl/@@download/$ctrl.bam"
done
samtools index -M *.bam

YAML Configuration¶

Save the following as demoHistoneChIPSeq.yaml in the directory containing the BAM files:

experimentName: demoHistoneChIPSeq

genomeParams:
  name: hg38
  chromosomes: [chr21, chr22]
  excludeForNorm: [chrX, chrY]

inputParams:
  samples:
    - name: H3K27ac_1
      path: ENCFF793ZHL.bam
    - name: H3K27ac_2
      path: ENCFF647VPO.bam
    - name: H3K27ac_3
      path: ENCFF809VKT.bam
    - name: H3K27ac_4
      path: ENCFF295EFL.bam
    - name: input_1
      path: ENCFF444WVG.bam
      role: control
    - name: input_2
      path: ENCFF619NYP.bam
      role: control
    - name: input_3
      path: ENCFF898LKJ.bam
      role: control
    - name: input_4
      path: ENCFF490MWV.bam
      role: control

Run Consenrich¶

% consenrich --config demoHistoneChIPSeq.yaml --verbose

The run will output a state bigWig, uncertainty bigWig, ROCCO narrowPeak file, and ROCCO gappedPeak file. With the current package version, the principal output files are:

demoHistoneChIPSeq_consenrich_state.VERSION.bw
demoHistoneChIPSeq_consenrich_uncertainty.VERSION.bw
consenrichOutput_demoHistoneChIPSeq_state.VERSION_rocco.narrowPeak
consenrichOutput_demoHistoneChIPSeq_state.VERSION_rocco.gappedPeak

The tags switch from state to stateShrunk and from uncertainty to uncertaintyShrunk for outputParams.stateShrinkageEnabled. To retain human-readable bedGraph files, set outputParams.bedGraphEnabled: true in the configuration.

ATAC-seq Demo¶

This demo estimates a consensus chromatin accessibility signal, uncertainty, and peaks from ten ENCODE ATAC-seq alignments.

Download Alignments¶

encodeFiles=https://www.encodeproject.org/files
for file in ENCFF326QXM ENCFF447ZRG ENCFF462RHM ENCFF495DQP ENCFF497QOS ENCFF632MBC ENCFF687QML ENCFF767FGV ENCFF919PWF ENCFF949CVL; do
    wget "$encodeFiles/$file/@@download/$file.bam"
done
samtools index -M *.bam

YAML Configuration¶

Save the following as atacDemo.yaml:

experimentName: atacDemo

genomeParams:
  name: hg38
  # chromosomes: [chr11, chr22] # faster: Uncomment and specify chromosomes
  excludeForNorm: [chrX, chrY]

inputParams:
  bamFiles:
    - ENCFF326QXM.bam
    - ENCFF447ZRG.bam
    - ENCFF462RHM.bam
    - ENCFF495DQP.bam
    - ENCFF497QOS.bam
    - ENCFF632MBC.bam
    - ENCFF687QML.bam
    - ENCFF767FGV.bam
    - ENCFF919PWF.bam
    - ENCFF949CVL.bam

Run Consenrich¶

% consenrich --config atacDemo.yaml --verbose

Results¶

Broad Mark ChIP-seq Demo¶

This demo estimates a consensus H3K4me1 signal from six ENCODE ChIP-seq experiments with matched input controls.

Download Alignments¶

encodeFiles=https://www.encodeproject.org/files

while read -r accession output; do
    wget -O "$output" "$encodeFiles/$accession/@@download/$accession.bam"
done <<'EOF'
ENCFF451FGF ENCFF451FGF_ENCSR724MJX_treatment_GRCh38.bam
ENCFF660DDB ENCFF660DDB_ENCSR438QZN_treatment_GRCh38.bam
ENCFF828GWI ENCFF828GWI_ENCSR564QBS_treatment_GRCh38.bam
ENCFF671AAF ENCFF671AAF_ENCSR299NYB_treatment_GRCh38.bam
ENCFF366HMH ENCFF366HMH_ENCSR817JNE_treatment_GRCh38.bam
ENCFF392MXC ENCFF392MXC_ENCSR485LPA_treatment_GRCh38.bam
ENCFF013ION ENCFF013ION_ENCSR109IWL_control_GRCh38.bam
ENCFF730PAF ENCFF730PAF_ENCSR526EXI_control_GRCh38.bam
ENCFF971XWL ENCFF971XWL_ENCSR945LPX_control_GRCh38.bam
ENCFF271VZL ENCFF271VZL_ENCSR261PLD_control_GRCh38.bam
ENCFF944AYC ENCFF944AYC_ENCSR061EXX_control_GRCh38.bam
ENCFF273EYI ENCFF273EYI_ENCSR120ITZ_control_GRCh38.bam
EOF

samtools index -M *.bam

YAML Configuration¶

Save the following as bigH3K4me1Demo.yaml:

experimentName: bigH3K4me1Demo

genomeParams:
  name: hg38
  chromosomes: [chr21, chr22]
  excludeForNorm: [chrX, chrY]

inputParams:
  samples:
    - name: ENCSR724MJX_H3K4me1
      path: ENCFF451FGF_ENCSR724MJX_treatment_GRCh38.bam
    - name: ENCSR438QZN_H3K4me1
      path: ENCFF660DDB_ENCSR438QZN_treatment_GRCh38.bam
    - name: ENCSR564QBS_H3K4me1
      path: ENCFF828GWI_ENCSR564QBS_treatment_GRCh38.bam
    - name: ENCSR299NYB_H3K4me1
      path: ENCFF671AAF_ENCSR299NYB_treatment_GRCh38.bam
    - name: ENCSR817JNE_H3K4me1
      path: ENCFF366HMH_ENCSR817JNE_treatment_GRCh38.bam
    - name: ENCSR485LPA_H3K4me1
      path: ENCFF392MXC_ENCSR485LPA_treatment_GRCh38.bam

    - name: ENCSR109IWL_input_for_ENCSR724MJX
      path: ENCFF013ION_ENCSR109IWL_control_GRCh38.bam
      role: control
    - name: ENCSR526EXI_input_for_ENCSR438QZN
      path: ENCFF730PAF_ENCSR526EXI_control_GRCh38.bam
      role: control
    - name: ENCSR945LPX_input_for_ENCSR564QBS
      path: ENCFF971XWL_ENCSR945LPX_control_GRCh38.bam
      role: control
    - name: ENCSR261PLD_input_for_ENCSR299NYB
      path: ENCFF271VZL_ENCSR261PLD_control_GRCh38.bam
      role: control
    - name: ENCSR061EXX_input_for_ENCSR817JNE
      path: ENCFF944AYC_ENCSR061EXX_control_GRCh38.bam
      role: control
    - name: ENCSR120ITZ_input_for_ENCSR485LPA
      path: ENCFF273EYI_ENCSR120ITZ_control_GRCh38.bam
      role: control

Run Consenrich¶

% consenrich --config bigH3K4me1Demo.yaml --verbose

Results¶

Consenrich for Single-Cell¶

In this demo, we use publicly available scATAC-seq data to infer signal and peaks for Goblet, Early Colonocyte, Colonocyte, and BEST4 cells in the ascending colon. Note that Consenrich can also accept data from other single-cell epigenomic assays, e.g., scCut&Tag.

Paper: Single-cell analysis of chromatin accessibility in the human intestine identifies regulatory programs and clarifies genetic associations in Crohn’s disease
GEO accession: GSE273194
Cell annotations: Zenodo

Donor	GEO sample	Goblet	Early Colonocyte	Colonocyte	BEST4
P3	GSM8423799	412	690	1277	46
P4	GSM8423800	647	983	1125	56
P6	GSM8423801	893	1355	1183	54
P12	GSM8423805	336	591	811	15
P13	GSM8423806	403	827	1034	74

Download Fragments and Cell Labels¶

This will be about 6 GB of data. You can remove donors and adapt the YAML config below as needed to reduce size.

In your terminal (assuming bash):

geoFiles=https://ftp.ncbi.nlm.nih.gov/geo/samples/GSM8423nnn

# you can use wget or curl, depending
while read -r accession donor; do
    fragmentFile="${accession}_${donor}-AC-CDnonInf_fragments.tsv.gz"
    curl -fLO "$geoFiles/$accession/suppl/$fragmentFile"
    curl -fLO "$geoFiles/$accession/suppl/$fragmentFile.tbi.gz"
    gzip -dc "$fragmentFile.tbi.gz" > "$fragmentFile.tbi"
done <<'EOF'
GSM8423799 P3
GSM8423800 P4
GSM8423801 P6
GSM8423805 P12
GSM8423806 P13
EOF

curl -fL -o scATAC_cell_metadata.csv \
    'https://zenodo.org/records/21508422/files/scATAC_cell_metadata.csv?download=1'

YAML Configuration¶

Select an input format to view its configuration. Each option uses five donors, four cell groups, and CPM normalization. The Pott downloads above match the first option.

See the single-cell API reference for field types and defaults.

Fragments + table

Use fragment files and the downloaded CSV cell labels.

scATACFragments.yaml¶

experimentName: scATACDemo
genomeParams:
  name: hg38
  chromosomes: [chr21, chr22]
countingParams:
  normMethod: CPM
scParams:
  defaultCountMode: center
  cellMetadata:
    path: scATAC_cell_metadata.csv
    barcodeColumn: barcode
    sampleColumn: Sample_ID
    groupColumn: cluster_label1
  fits:
    - {name: Goblet, selectGroups: [AC_Goblet]}
    - {name: EarlyColonocyte, selectGroups: [Early_Colonocyte]}
    - {name: Colonocyte, selectGroups: [Colonocyte]}
    - {name: BEST4, selectGroups: [BEST4]}
inputParams:
  samples:
    - name: P3
      path: GSM8423799_P3-AC-CDnonInf_fragments.tsv.gz
      metadataSample: P3-AC-CDnonInf
    - name: P4
      path: GSM8423800_P4-AC-CDnonInf_fragments.tsv.gz
      metadataSample: P4-AC-CDnonInf
    - name: P6
      path: GSM8423801_P6-AC-CDnonInf_fragments.tsv.gz
      metadataSample: P6-AC-CDnonInf
    - name: P12
      path: GSM8423805_P12-AC-CDnonInf_fragments.tsv.gz
      metadataSample: P12-AC-CDnonInf
    - name: P13
      path: GSM8423806_P13-AC-CDnonInf_fragments.tsv.gz
      metadataSample: P13-AC-CDnonInf

% consenrich --config scATACFragments.yaml --verbose2

Fragments + H5AD

Use fragment files and labels stored in scATAC_cell_metadata.h5ad. In this case, the label file supplies .obs only.

For an H5AD option, install the optional dependency:

python -m pip install "consenrich[singlecell]"

See the single-cell API reference for H5AD input requirements and metadata fields.

scATACFragmentsH5ADMetadata.yaml¶

experimentName: scATACDemo
genomeParams:
  name: hg38
  chromosomes: [chr21, chr22]
countingParams:
  normMethod: CPM
scParams:
  defaultCountMode: center
  cellMetadata:
    path: scATAC_cell_metadata.h5ad
    barcodeColumn: barcode
    sampleColumn: Sample_ID
    groupColumn: cluster_label1
  fits:
    - {name: Goblet, selectGroups: [AC_Goblet]}
    - {name: EarlyColonocyte, selectGroups: [Early_Colonocyte]}
    - {name: Colonocyte, selectGroups: [Colonocyte]}
    - {name: BEST4, selectGroups: [BEST4]}
inputParams:
  samples:
    - name: P3
      path: GSM8423799_P3-AC-CDnonInf_fragments.tsv.gz
      metadataSample: P3-AC-CDnonInf
    - name: P4
      path: GSM8423800_P4-AC-CDnonInf_fragments.tsv.gz
      metadataSample: P4-AC-CDnonInf
    - name: P6
      path: GSM8423801_P6-AC-CDnonInf_fragments.tsv.gz
      metadataSample: P6-AC-CDnonInf
    - name: P12
      path: GSM8423805_P12-AC-CDnonInf_fragments.tsv.gz
      metadataSample: P12-AC-CDnonInf
    - name: P13
      path: GSM8423806_P13-AC-CDnonInf_fragments.tsv.gz
      metadataSample: P13-AC-CDnonInf

% consenrich --config scATACFragmentsH5ADMetadata.yaml --verbose2

H5AD cells

This case uses one raw cell-by-bin count matrix per donor, such as P3.cells.h5ad. matrixRows: cells sums selected cells per donor and cell group into internal pseudobulks.

For an H5AD option, install the optional dependency:

python -m pip install "consenrich[singlecell]"

See the single-cell API reference for H5AD input requirements and metadata fields.

scATACCells.yaml¶

experimentName: scATACDemo
genomeParams:
  name: hg38
  chromosomes: [chr21, chr22]
countingParams:
  normMethod: CPM
scParams:
  fits:
    - {name: Goblet, selectGroups: [AC_Goblet]}
    - {name: EarlyColonocyte, selectGroups: [Early_Colonocyte]}
    - {name: Colonocyte, selectGroups: [Colonocyte]}
    - {name: BEST4, selectGroups: [BEST4]}
inputParams:
  samples:
    - name: P3
      path: P3.cells.h5ad
      matrixRows: cells
      cellMetadata: {path: P3.cells.h5ad, groupColumn: cluster_label1}
    - name: P4
      path: P4.cells.h5ad
      matrixRows: cells
      cellMetadata: {path: P4.cells.h5ad, groupColumn: cluster_label1}
    - name: P6
      path: P6.cells.h5ad
      matrixRows: cells
      cellMetadata: {path: P6.cells.h5ad, groupColumn: cluster_label1}
    - name: P12
      path: P12.cells.h5ad
      matrixRows: cells
      cellMetadata: {path: P12.cells.h5ad, groupColumn: cluster_label1}
    - name: P13
      path: P13.cells.h5ad
      matrixRows: cells
      cellMetadata: {path: P13.cells.h5ad, groupColumn: cluster_label1}

% consenrich --config scATACCells.yaml --verbose2

H5AD pseudobulks

This case uses one raw pseudobulk-by-bin count matrix per donor, such as P3.pseudobulks.h5ad. matrixRows: pseudobulks selects one prepared pseudobulk per donor and cell group.

For an H5AD option, install the optional dependency:

python -m pip install "consenrich[singlecell]"

See the single-cell API reference for H5AD input requirements and metadata fields.

scATACPseudobulks.yaml¶

experimentName: scATACDemo
genomeParams:
  name: hg38
  chromosomes: [chr21, chr22]
countingParams:
  normMethod: CPM
scParams:
  fits:
    - {name: Goblet, selectGroups: [AC_Goblet]}
    - {name: EarlyColonocyte, selectGroups: [Early_Colonocyte]}
    - {name: Colonocyte, selectGroups: [Colonocyte]}
    - {name: BEST4, selectGroups: [BEST4]}
inputParams:
  samples:
    - name: P3
      path: P3.pseudobulks.h5ad
      matrixRows: pseudobulks
      cellMetadata: {path: P3.pseudobulks.h5ad, groupColumn: cluster_label1}
    - name: P4
      path: P4.pseudobulks.h5ad
      matrixRows: pseudobulks
      cellMetadata: {path: P4.pseudobulks.h5ad, groupColumn: cluster_label1}
    - name: P6
      path: P6.pseudobulks.h5ad
      matrixRows: pseudobulks
      cellMetadata: {path: P6.pseudobulks.h5ad, groupColumn: cluster_label1}
    - name: P12
      path: P12.pseudobulks.h5ad
      matrixRows: pseudobulks
      cellMetadata: {path: P12.pseudobulks.h5ad, groupColumn: cluster_label1}
    - name: P13
      path: P13.pseudobulks.h5ad
      matrixRows: pseudobulks
      cellMetadata: {path: P13.pseudobulks.h5ad, groupColumn: cluster_label1}

% consenrich --config scATACPseudobulks.yaml --verbose2

Run the command shown in the selected tab from the directory containing its YAML file and inputs. Outputs are written under scATACDemo/<cell group>/.

Results¶

Browser snapshots for each cell type (uncertainty tracks omitted for spacing). Signac and ArchR are included here as familiar references, both run with generic defaults.

Quickstart + Usage¶

Getting Started: H3K27ac ChIP-seq¶

Input Data¶

Download Alignments¶

YAML Configuration¶

Run Consenrich¶

ATAC-seq Demo¶

Download Alignments¶

YAML Configuration¶

Run Consenrich¶

Results¶

Broad Mark ChIP-seq Demo¶

Download Alignments¶

YAML Configuration¶

Run Consenrich¶

Results¶

Consenrich for Single-Cell¶

Download Fragments and Cell Labels¶

YAML Configuration¶

Results¶

Goblet¶

Early Colonocyte¶

Colonocyte¶

BEST4¶