Skip to content

daklab/tealeaf

Folders and files

NameName
Last commit message
Last commit date

Latest commit

 

History

98 Commits
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 

Repository files navigation

tealeaf ( previously LeafCutterITI)

By Xingpei Zhang and David A. Knowles

Citations:

Alamancos, G. P., Pagès, A., Trincado, J. L., Bellora, N., & Eyras, E. (2015). Leveraging transcript quantification for fast computation of alternative splicing profiles. RNA , 21(9), 1521–1531. https://doi.org/10.1261/rna.051557.115

Garrido-Martín, D., Palumbo, E., Guigó, R., & Breschi, A. (2018). ggsashimi: Sashimi plot revised for browser-and annotation-independent splicing visualization. PLoS computational biology, 14(8), e1006360.

Li, Y. I., Knowles, D. A., Humphrey, J., Barbeira, A. N., Dickinson, S. P., Im, H. K., & Pritchard, J. K. (2018). Annotation-free quantification of RNA splicing using LeafCutter. Nature Genetics, 50(1), 151–158. https://doi.org/10.1038/s41588-017-0004-9

Patro, R., Duggal, G., Love, M. I., Irizarry, R. A., & Kingsford, C. (2017). Salmon provides fast and bias-aware quantification of transcript expression. Nature Methods, 14(4), 417–419. https://doi.org/10.1038/nmeth.4197

Requirements (versions used for development)

  • python (v3.10.11)

Python Dependencies

  • numpy
  • pandas
  • pyranges
  • scipy
  • scikit-learn
  • pytorch
  • tqdm
  • pyro-ppl
  • matplotlib

R Dependencies (for tealeaf-ggsashimi)

  • R
  • ggplot2
  • data.table
  • gridExtra

Additional Requirement for isoform quantification

  • salmon (v1.10.0)

Other dependencies for Leafcutter as listed in https://github.com/davidaknowles/leafcutter/tree/master, especially for Leafcutter_ds

tealeaf

A modified version of Leafcutter that detects and analyzes alternative splicing events by quantifying excised introns by utilizing isoform abundance and transcriptome annotation. Can also be install as a command line tool with pip install tealeaf

tealeaf_Workflow

There are five main command-line tools in tealeaf:

  • tealeaf_map_gen
  • tealeaf_clustering (for bulk & bulk-like single-cell or pseudobulk)
  • tealeaf_sc (for single-cell)
  • tealeaf-ggsashimi (utlize adapted ggsashimi for sashimi plotting)
  • tealeaf-test (Dirichlet-multinomial differential splicing test with empirical-null calibration)

tealeaf_map_gen

usage: python tealeaf_map_gen.py [-a/--annot] [--annot_source] [-o/--outprefix] 
                     [--maxintronlen] [--minintronlen] [-v/--virtual_intron] [--single_cell]
or when install with pip
tealeaf-map [-a/--annot] [--annot_source] [-o/--outprefix] [--maxintronlen]
                      [--minintronlen] [-v/--virtual_intron] [--single_cell]


Mandatory parameters:

-a, --annot     The transcriptome annotation gtf file for tealeaf-map to run with 

Optional Parameters:

--annot_source          The annotation source for the annotation, currently support Gencode and Stringtie 
                        (default: gencode)

-o, --outprefix         The prefix for output files (default: Leafcutter_)

--maxintronlen          The maximum allowed intron length for introns (default: 5,000,000)

--minintronlen          The minimum allowed intron length for introns (default: 50)

--no_quality_control  A flag on whether not to remove pseudogene, and decay transcript

-v, --virtual_intron    A flag on whether to compute virtual intron that can be used to capture
                        AFE and ALE usage, a testing feature

--single_cell           Whether to build matrices for isoform to intron and exon, required if dealing with\
                        single cell data from alevin-fry (default: True)

tealeaf_clustering

usage: python tealeaf_clustering.py [--map] [--count_files] [--connect_file] [-a/--annot]
                    [--cluster_def] [-o/--outprefix] [--use_TPM] [--samplecutoff]
                    [--introncutoff] [-m/--minclucounts] [-r/--mincluratio] [--normalization_scale]
                    [--read_len] [--not_paired_end] [--overhang] [--sizing_factor]
or when install with pip
tealeaf-cluster [--map] [--count_files] [--connect_file] [-a/--annot]
                    [--cluster_def] [-o/--outprefix] [-n/--normalization] [--samplecutoff]
                    [--introncutoff] [-m/--minclucounts] [-r/--mincluratio] [--normalization_scale]
                    [--read_len] [--not_paired_end] [--overhang] [--sizing_factor]


Mandatory parameters:
--map             The isoforms to introns map generated from tealeaf-map  

--count_files     A txt file that contain the sample names 

--connect_file    The intron-exon connectivity file generated from tealeaf-map 

-a, --annot       The transcriptome annotation gtf file 


Optional Parameters:
--cluster_def           The definition used for cluster refinement, three def available, 1: overlap, 2: overlap+share_intron_splice_site, 
                        3: overlap+share_intron_splice_site+shared_exon_splice_site (default: 3)

-o, --outprefix         The prefix for output files (default: Leafcutter_)

--use_TPM               A flag on whether to use TPM or normalized count

--preprocessed          A flag on Whether the files provided are already normalized, mainly for rerunning the pipeline and don't 
                        perform normalization again
--normalization_scale   The mode use for normaliztion, whether the count/TPM scale is based on junction count simulation, local (gene level) or global level,
                        can only input junction, local, or global (default: junction)

--samplecutoff          Minimum Normalized count/TPM for an intron in a sample to count as exist (default: 0)

--introncutoff          Minimum Normalized count/TPM for an intron to count as exist(default 5)

--m, --minclucounts     Minimum Normalized count/TPM to support a cluster (default: 30)

-r, --mincluratio       Minimum fraction of reads in a cluster that supports an intron (default: 0.01)

--read_len              The read length of sequencing data, use to simulate junction count, only work when normalization_scale = "junction"

--not_paired_end        Whether the reads are not paired end use to simulate junction count, only work when normalization_scale = "junction"
                        (default: False)
--overhang              The oeverhand that would like to use, could be set to zero, use to simulate junction count, 
                        only work when normalization_scale = junction (default: 2)

--sizing_factor         The sizing factor for junction simulation normalization to better calibrate the p-values (default: 1)

tealeaf_sc

usage: python tealeaf_sc.py [--alevin_dir] [--salmon_ref] [--ref_dir] [--barcodes_clusters] [--pseudobulk_samples]
                          [-n/--num_cell] [-k/--num_bootstrapping] [--min_eq] [--group_method] [--ref_prefix]
                          [--thread] [--cluster_def] [-o/--outprefix] [-n/--normalization] [--samplecutoff] 
                          [--introncutoff] [-m/--minclucounts] [-r/--mincluratio] [--preprocessed]
                          [--read_len] [--not_paired_end] [--overhang] [--sizing_factor]


or when install with pip
tealeaf-sc [--alevin_dir] [--salmon_ref] [--ref_dir] [--barcodes_clusters] [--pseudobulk_samples]
                          [-n/--num_cell] [-k/--num_bootstrapping] [--min_eq] [--group_method] [--ref_prefix]
                          [--thread] [--cluster_def] [-o/--outprefix] [--normalization_scale] [--samplecutoff] 
                          [--introncutoff] [-m/--minclucounts] [-r/--mincluratio] [--preprocessed]
                          [--read_len] [--not_paired_end] [--overhang] [--sizing_factor]


Mandatory parameters:

--alevin_dir            The directory for alevin results, the file should contain the eq matrix and other files

--salmon_ref            The reference used for salmon index, The salmon reference,  maybe spliceu or splicei

--ref_dir               Tealeaf reference directory, which should contain the matrices for isoform to intron and exon

--barcodes_clusters     The file that records which barcodes belong to which cluster/cell type in the format 'barcode,cluster' 
                        this file will be used to generate pseudobulk samples 

--pseudobulk_samples    A txt file with barcodes to pseudobulk sample are expected in format 'barcode pseudobulk_ample', if \
                        this option != None, then it will overwrite the input to --barcodes_clusters, and use the file in this option \
                        for computation. Only one of barcodes_clusters or pseudobulk_samples is required

Optional Parameters:
--ref_prefix            The prefix that is used to generate isoform to intron map using
                        tealeaf-map (default: '')

--n,--num_cell          The number of cell/barcode that you would like to include in a pseudobulk sample, cluster/cell type that has fewer
                        cell/barcodes than this number will not included in the computation (default: 100)

-k,--num_bootstrapping  The number of bootstrapping samples generated for each cluster/cell type if using bootstrapping to generate pseudobulk sample (default: 30)

--min_eq                Minimum count for each eq class for it to be included in the EM (default: 5)

--pseudobulk_method     The pseudobulk sample generate method, could be metacells or bootstrapping (default: metacells)

--cluster_def           The definition used for cluster refinement, three def available, 1: overlap, 2: overlap+share_intron_splice_site, 
                        3: overlap+share_intron_splice_site+shared_exon_splice_site (default: 3)

-o, --outprefix         The prefix for output files (default: leafcutter_)

--thread                The number of threads used for parallel computation, should not be too large to avoid crash (default: 8)

--use_TPM               A flag on whether to use TPM or normalized count

--preprocessed          A flag on whether pseudobulk generation and EM were done, if true, then the pipeline starts from counting intron (default: False)

-v,--with_virtual       A flag on whether the map that contain virtual intron to capture AFE and ALE

--samplecutoff          Minimum Normalized count/TPM for an isoform in a sample to count as exist (default: 0.1)

--introncutoff          Minimum Normalized count/TPM for an intron to count as exist(default: 80)

--m, --minclucounts     Minimum Normalized count/TPM to support a cluster (default: 100)

-r, --mincluratio       Minimum fraction of reads in a cluster that supports an intron (default 0.01)

--normalization_scale   The mode use for normaliztion, whether the count/TPM scale is based on junction count simulation, local (gene level) or global level,
                        can only input junction, local, or global (default: junction)

--read_len              The read length of sequencing data, use to simulate junction count, only work when normalization_scale = "junction" (default: 100)

--not_paired_end        Whether the reads are not paired end use to simulate junction count, only work when normalization_scale = "junction" (default: False)
--overhang              The oeverhand that would like to use, could be set to zero, use to simulate junction count, 
                        only work when normalization_scale = junction (default: 2)

--sizing_factor         The sizing factor for junction simulation normalization to better calibrate the p-values (default: 1)

tealeaf_test

usage: tealeaf-test [-h] [-0 BASELINE_GROUP] [-o OUTPUT_PREFIX]
                    [-s MAX_CLUSTER_SIZE] [-i MIN_SAMPLES_PER_INTRON]
                    [-g MIN_SAMPLES_PER_GROUP] [-c MIN_COVERAGE]
                    [-u MIN_UNIQUE_VALS] [--init INIT] [-p NUM_THREADS]
                    [--lbfgs_max_iter LBFGS_MAX_ITER]
                    [--lbfgs_history_size LBFGS_HISTORY_SIZE]
                    [--pseudobulk_size PSEUDOBULK_SIZE] [--no_pseudobulk]
                    [--observed_reps OBSERVED_REPS] [--null_reps NULL_REPS]
                    [--seed SEED]
                    counts_file [groups_file]

tealeaf-test runs a Dirichlet-multinomial differential splicing test on an intron-by-sample count matrix. It can also run an empirical-null calibration step by repeatedly permuting group labels and using the resulting null log-likelihood ratios to calibrate the observed test statistics.

The Dirichlet-multinomial statistical model and core test implementation are based on the LeafCutter2 stat-test code from the Knowles lab LeafCutter repository: https://github.com/davidaknowles/leafcutter/tree/leafcutter2/leafcutter/differential_splicing. Tealeaf adds package integration and empirical-null calibration around this test.

Input requirements:

  • counts_file: whitespace-delimited intron-by-sample count matrix. Row identifiers must use chr:start:end:cluster format. Columns are samples or already-created pseudobulk samples.
  • groups_file: optional whitespace-delimited sample metadata file with at least two columns: sample name and group label. Extra columns are treated as confounders. If no group file is provided, group labels are inferred from sample names such as astrocyte_0.bam or astrocyte_0_quant_normalized.sf.

Common options:

  • --no_pseudobulk: use the existing count columns directly. This is the recommended option when the input has already been pseudobulked or already represents sample-level observations.
  • --pseudobulk_size: create pseudobulk replicates before testing. Leave this at 0 when using precomputed pseudobulk files.
  • --observed_reps: number of observed pseudobulk replicates to test. For direct sample mode, this is usually 1.
  • --null_reps: number of empirical-null label permutations. Use a larger value for final analyses and a small value for quick testing.
  • -s, --max_cluster_size: skip clusters with more introns than this threshold.
  • -i, --min_samples_per_intron: require each retained intron to be used in at least this many samples.
  • -g, --min_samples_per_group: require at least this many covered samples in at least two groups.
  • -c, --min_coverage: minimum cluster coverage used for valid-sample checks.
  • -p, --num_threads: number of worker processes.

The main cluster output includes the model p-value (p), BH-adjusted model p-value (p.adjust), empirical-null calibrated p-value (empirical_null_p), and BH-adjusted empirical-null p-value (empirical_null_p.adjust). The empirical-null adjusted p-value is computed by applying BH correction directly to empirical_null_p.

Example:

tealeaf-test muris_refined_cluster sample_groups.tsv \
  -o muris_tealeaf_test --no_pseudobulk --null_reps 2 -p 4

For a final run, increase --null_reps to the desired number of permutations:

tealeaf-test muris_refined_cluster sample_groups.tsv \
  -o muris_tealeaf_test --no_pseudobulk --null_reps 100 -p 8

Main outputs:

  • {output_prefix}_cluster_significance.tsv: one row per cluster, including raw, BH-adjusted, empirical-null, and BH-adjusted empirical-null p-values. Failed or skipped clusters are retained with a status reason.
  • {output_prefix}_observed_cluster_all.tsv: observed test results for each observed replicate.
  • {output_prefix}_empirical_null_cluster_all.tsv: per-cluster results from empirical-null permutations.
  • {output_prefix}_empirical_null_loglr.tsv: empirical-null log-likelihood-ratio values used for calibration.

Clusters are considered testable when they pass the coverage, intron-use, group-size, and cluster-size filters and the model fit succeeds. Non-testable clusters are kept in the output with the reason in the status column.

tealeaf_ggsashimi

usage: tealeaf-ggsashimi [--intron INTRON] [--exon EXON]
                         [--strand_info STRAND_INFO]
                         [--data_type DATA_TYPE] [--normalized NORMALIZED]
                         [--aggregation AGGREGATION]
                         -c COORDINATES -g GTF
                         [-o OUT_PREFIX] [-S OUT_STRAND]
                         [-M MIN_COVERAGE] [-j JUNCTIONS_BED]
                         [-s STRAND] [--shrink] [-A AGGR]
                         [--alpha ALPHA] [-P PALETTE]
                         [--fix-y-scale] [--height HEIGHT]
                         [--ann-height ANN_HEIGHT] [--width WIDTH]
                         [--base-size BASE_SIZE]
                         [-F OUT_FORMAT] [-R OUT_RESOLUTION]
                         [--debug-info]

tealeaf-ggsashimi creates sashimi plots from Tealeaf intron and exon count tables for a selected genomic region. It is adapted from ggsashimi and uses Tealeaf count outputs directly instead of BAM files.

Required inputs for a normal Tealeaf plot:

  • --intron: Tealeaf intron count table.
  • --exon: Tealeaf exon count table.
  • --strand_info: intron-exon connectivity file from tealeaf_map_gen; used to assign strand information to introns and nearby exons.
  • -c, --coordinates: genomic interval to plot, in chr:start-end format.
  • -g, --gtf: annotation GTF file. Exon annotation is sufficient.

Input interpretation options:

  • --data_type: input type used to parse sample names. Use bulk or sc [default: bulk].
  • --normalized: whether input column names follow normalized Tealeaf naming [default: True].
  • --aggregation: whether to aggregate Tealeaf values by cell type/sample group before plotting [default: True]. Use either this option or -A/--aggr depending on whether aggregation should happen before or during ggsashimi-style overlay plotting.

Plot and output options:

  • -o, --out-prefix: output plot prefix [default: sashimi].
  • -F, --out-format: output format: pdf, svg, png, jpeg, or tiff [default: pdf].
  • -R, --out-resolution: raster output resolution in PPI [default: 300].
  • -M, --min-coverage: minimum junction count required before drawing a junction arc [default: 1].
  • -s, --strand: strand specificity: NONE, SENSE, ANTISENSE, MATE1_SENSE, or MATE2_SENSE [default: NONE].
  • -S, --out-strand: for stranded plots, choose both, plus, or minus [default: both].
  • --shrink: shrink intronic regions for a more compact display.
  • -A, --aggr: overlay aggregation function: mean, median, mean_j, or median_j.
  • --alpha: density histogram transparency [default: 0.5].
  • -P, --palette: color palette file; first column should contain R color names or hex colors.
  • --fix-y-scale: use the same y-axis scale across signal tracks.
  • --height, --ann-height, --width, --base-size: control signal-track height, annotation height, plot width, and base font size.
  • -j, --junctions-bed: optionally write the plotted junctions to a BED file.
  • --debug-info: print system/debug information useful for troubleshooting R and plotting dependencies.

Example:

tealeaf-ggsashimi \
  --intron sample_intron_counts.tsv \
  --exon sample_exon_counts.tsv \
  --strand_info sample_intron_exon_connectivity.tsv \
  -c chr1:100000-120000 \
  -g annotation.gtf \
  -o sample_sashimi \
  -F pdf --data_type bulk --normalized True --aggregation True

Main output:

  • {out_prefix}.pdf when --strand NONE and -F pdf are used.
  • {out_prefix}_plus.{format} and/or {out_prefix}_minus.{format} for stranded output when --strand is not NONE.

Detailed Tutorial to run the tealeaf

In this tutorial, we walk through all the steps to run the tealeaf pipeline. For each step, we discuss the possible parameters that can be changed, how to do so and the considerations involved in each of the parameters. Finally, we show example inputs and outputs of each step (with column explanations) so the user knows what to expect and can make custom files as needed.

Step 0: Transcriptome annotation download or generation and Salmon isoform quantification

Example human transcriptome annotation can be downloaded from https://www.gencodegenes.org/human/

Step 1: Isoform to intron map generation

In this step, tealeaf_map_gen will be used to generate a map that contains information about which isoform is generated by splicing which introns. The map will also contain information about which exon is in which isoform. This step only needs to run once for each unique transcriptome annotation gtf file.

Sample run:

python tealeaf_map_gen -a gencode.v45.annotation.gtf --annot_source gencode -o sample_run_ --maxintronlen 5000000 --minintronlen 50 -v False          

Depending on the setting, two or four files will be generated.

  • {out_prefix}isoform_intron_map.tsv
  • {out_prefix}intron_exon_connectivity.tsv
  • {out_prefix}isoform_intron_map_with_virtual.tsv
  • {out_prefix}intron_exon_connectivity_with_virtual.tsv

where with_virtual mean virtual intron was used to capture all annotated AFE and ALE (a testing feature).

if annotation_source='gencode', an additional file will be generated to give out information about the possible isoform type that can be generated by splicing out each intron

  • {out_prefix}intron_source_map.tsv (also a testing feature)

A record file that contains the parameters will also be generated

When --single_cell == True, five additional files will be generated. Two for sparse matrices in npz format, rows are isoforms, and columns are introns or exons. Three txt files record the row and column names. These file is essential for tealeaf-sc.

Step 2: Salmon isoform quantification

tealeaf utlized pseudoalignment method Salmon for bulk and preprocessed pseudobulk data. For usage of Salmon please refer https://salmon.readthedocs.io/en/latest/salmon.html

For single-cell data, tealeaf uses the alevin-fry pipeline from Salmon. The usage of alevin-fry please refer https://alevin-fry.readthedocs.io/en/latest/. Specific notices, please use -d, --dump-eqclasses flag when using alevin-fry quant to obtain the eqclass matrix. Also, t2t mapping should be used instead of normal t2g mapping. t2t mapping can be easily obtained by replacing the gene col in t2g file with the transcripts.

In the rest of the tutorial, we assume RNA-seq data aligned to the transcriptome using Salmon or Alevin-fry.

Step 2.1: Single-cell clustering after alevin-fry

For single-cell data, after pseudoaligment, we will need to process the data and obtain a barcodes to clusters/cell_types csv file that have row in format 'barcode,cluster/cell_type'. There are different single-cell analysis tool can achieve this goal. For examples, Seurat and Scanpy. Any analysis tool could work as long as the barcodes to clusters/cell_types csv file is provided. For our analysis, we used Scanpy and tutorial for cell clustering with Scanpy could be found at https://scanpy.readthedocs.io/en/stable/tutorials/basics/clustering.html. After the clustering and cell labeling, the barcodes to clusters/cell_types could be export like adata.obs[['cell_barcodes', 'cluster_name']].to_csv('barcode_to_cluster.csv', index = False, header = None)

Step 3.1: tealeaf clustering for bulk or pseudobulk data

For this step, we assume the data we are processing are bulk or pseudobulk data, we will need the file generated from step 1, the file path and the name for the isoform quantification files generated by Salmon, and the transcriptome annotation

Sample run:

python tealeaf_clustering --map transcript_intron_map.tsv --count_files quantification_files.txt --connect_file intron_exon_connectivity.tsv -a gencode.v45.annotation.gtf --cluster_def 3 \
                                --normalization True -o sample_run_ --minclucounts 30 --mincluratio 0.01

Two main output files will be obtained:

  • {outprefix}refined_cluster
  • {outprefix}ratio_count

sample {out_prefix}refined_cluster

sample1.sf sample2.sf sample3.sf sample4.sf sample5.sf sample6.sf
chr1:17055:17233:clu_1 21.1 13 18 20 17 12 
chr1:17055:17606:clu_1 4 11.4 12 7 2 0 5 
chr1:17368:17606:clu_1 127 132 128 55 93 90 68 
chr1:668593:668687:clu_2 3 11.3 1 3 4 4 8 
chr1:668593:672093:clu_2 11 16 23 2.5 3 20 9

These two files are equivalent to Leafcutter clustering numers.counts.gz and counts.gz. It is worth noticing that the normalized count or TPM is not necessarily an integer, but the normalized count will exhibit a count-like property.

Step 3.2: tealeaf clustering for single-cell data

For this step, we assume the data are single-cell data and that the outputs from Step 2 and Step 2.1 are available. tealeaf-sc creates pseudobulk samples from cells, estimates pseudobulk intron/exon usage, and then performs Tealeaf cluster refinement.

The required inputs are:

  • --alevin_dir: the alevin-fry output directory. It should contain files such as gene_eqclass.txt.gz, geqc_counts.mtx, and the other relevant alevin-fry quantification outputs.
  • --salmon_ref: the Salmon reference FASTA used to build the index, such as spliceu.fa or splicei.fa.
  • --ref_dir: the Tealeaf reference directory generated by tealeaf-map; this directory should contain the isoform-to-intron and isoform-to-exon matrices.
  • --barcodes_clusters: a barcode-to-cluster/cell-type file in barcode,cluster format. This is used to generate pseudobulk samples.
  • --pseudobulk_samples: optional barcode-to-pseudobulk file in barcode pseudobulk_sample format. If provided, this overrides --barcodes_clusters.

Common options:

  • -n, --num_cell: number of cells/barcodes per pseudobulk sample [default: 100].
  • -k, --num_bootstrapping: number of bootstrapping samples when using bootstrapping pseudobulk generation [default: 30].
  • --pseudobulk_method: pseudobulk method, either metacells or bootstrapping [default: metacells].
  • --cluster_def: cluster refinement definition. 1 uses overlap, 2 uses overlap plus shared intron splice site, and 3 uses overlap plus shared intron and exon splice sites [default: 3].
  • --normalization_scale: normalization scale, such as junction, snRNA, or global [default: junction].
  • --preprocessed: start from preprocessed pseudobulk/EM files and skip pseudobulk generation and EM.
  • --thread: number of threads used for computation [default: 8].

Sample run:

tealeaf-sc \
  --alevin_dir salmon/out_permit_know/quant_spliceu_t2t \
  --salmon_ref salmon_index/spliceu.fa \
  --ref_dir tealeaf_map_dir \
  --barcodes_clusters barcodes_clusters.txt \
  --cluster_def 3 \
  --normalization_scale junction \
  -o sample_run_ \
  --minclucounts 30 --mincluratio 0.01 -n 100

Similar to tealeaf-cluster, two main output files will be obtained:

  • {outprefix}refined_cluster
  • {outprefix}ratio_count

This step also generates supporting files such as barcodes_to_pseudobulk.txt, which records the pseudobulk assignments used during the computation.

Step 4: Differential splicing test with empirical-null calibration

The refined_cluster output from Step 3 can be tested with tealeaf-test. This command runs a Dirichlet-multinomial likelihood-ratio test for each cluster and can optionally calibrate the observed p-values using empirical-null label permutations.

The main input is an intron-by-sample count matrix where row IDs have the format chr:start:end:cluster. The Tealeaf {outprefix}refined_cluster file already follows this format. The optional group file is a whitespace-delimited table with at least two columns:

sample1 Control
sample2 Control
sample3 Treatment
sample4 Treatment

If no group file is supplied, tealeaf-test will try to infer group labels from the sample column names. For final analyses, an explicit group file is recommended.

For already pseudobulked or sample-level input, use --no_pseudobulk so the existing columns are used directly:

tealeaf-test sample_run_refined_cluster sample_groups.tsv \
  -o sample_run_tealeaf_test \
  --no_pseudobulk \
  --null_reps 100 \
  -p 8

For a quick test run, use fewer permutations:

tealeaf-test sample_run_refined_cluster sample_groups.tsv \
  -o sample_run_tealeaf_test_quick \
  --no_pseudobulk \
  --null_reps 2 \
  -p 2

Important filtering options:

  • -s, --max_cluster_size: skip clusters with more introns than this threshold.
  • -i, --min_samples_per_intron: require each retained intron to be observed in at least this many samples [default: 5].
  • -g, --min_samples_per_group: require at least this many valid samples in at least two groups [default: 3].
  • -c, --min_coverage: minimum cluster coverage used for valid-sample checks [default: 20].
  • --observed_reps: number of observed pseudobulk replicates to test [default: 1].
  • --null_reps: number of empirical-null label permutations.

Main outputs:

  • {output_prefix}_cluster_significance.tsv: cluster-level test table with p, p.adjust, empirical_null_p, empirical_null_p.adjust, and status.
  • {output_prefix}_observed_cluster_all.tsv: observed test results for each observed replicate.
  • {output_prefix}_empirical_null_cluster_all.tsv: per-cluster empirical-null permutation results.
  • {output_prefix}_empirical_null_loglr.tsv: empirical-null log-likelihood-ratio values used for calibration.

Clusters with status missing or successful are testable. Non-testable clusters remain in the output with a status message explaining which filter they failed.

About

No description, website, or topics provided.

Resources

License

Stars

0 stars

Watchers

1 watching

Forks

Releases

No releases published

Packages

 
 
 

Contributors

Languages