This repository contains the code to reproduce the results from the Birth Defects Biorespository (BDB) paper. Each subfolder of this repository represents a major analysis that was performed and will contain its own README.md with instructions on how to run the analysis along with the necessary code.
lrRNA-seq = long-read RNA-sequencing
Our paper uses hg38 as the human genome reference and GENCODE v44 primary assembly comprehensive gene annotations. Feel free to swap these out for reference files you are interested in using.
We assume that VCF files for cohort-level DNA sequencing data is present. There are numerous established pipelines for variant calling on DNA sequencing data. One example includes the GATK HaplotypeCaller pipeline, which we used to get our DNA variant calls. Note that our scripts assume that chromosome names in all VCF files begin with chr.
This analysis processes raw nanopore sequencing data into bam files. This analysis is required for all downstream analyses.
Please refer to the nanopore_data_processing folder for full details.
This analysis checks the agreement between variants called from lrRNA-seq and genome sequencing. This ensures that no sample swaps occured during processing of the samples.
Please refer to the dna_rna_check folder for full details.
This analysis discovers and quantifies full-length transcript isoforms using our lab's ESPRESSO tool. This analysis is required for all downstream analysis.
Please refer to the espresso folder for full details.
This analysis discovers and quantifies the usage of all splice junctions in a lrRNA-seq cohort. A large majority of the splicing analysis in this paper utilizes these splice junction usage estimates. See our paper methods for more details about how splice junction usage is quantified. This analysis is required for all downstream analysis.
Please refer to the sj_quant folder for full details.
This analysis preprocesses genome sequencing data (in our case, whole genome sequencing data) to prepare it for downstream analyses. This analysis is required for all downstream analysis.
Please refer to the preprocess_genome_sequencing folder for full details.
This analysis discovers common splicing variants and characterizes their molecular and phenotypic effects. Specifically, this analysis includes splice junction QTL (sjQTL) mapping, fine-mapping of sjQTL loci, sjQTL event classification, GWAS colocalization, CDK2AP/TLR1/CDK10 splicing event QTL analyses highlighted in the paper, and GTEx overlap analysis.
Please refer to the common_splicing_variant_analysis folder for full details.
This analysis retrieves the gnomAD allele frequency (AF) of each variant in the cohort's genome sequencing data. This will be used in the allele-specific and genetic disease diagnosis analyses.
Please refer to the get_gnomad_af folder for full details.
This analysis uses allele-specific splice junction analysis to discover, for each sample in the cohort, splice junctions that are differentially used between the two haplotypes of that sample. It also finds all proximal heterozygous variants near a differentially used splice junction that could explain the observed allele-specific usage of the splice junction. Additional downstream analysis include: comparison of allele specific splice junction associated variants to sjQTLs, SpliceAI prioritization of likely causal splicing mutations responsible for allele-specific splice junctions, analysis of SpliceAI prioritized variants across genes and the cohort, implication of SpliceAI prioritized variants in GWAS loci, and the CARD9 and PUF60 analyses shown in the paper.
Please refer to the allele_specific_analysis folder for full details.
This analysis integrates genetic variant data, lrRNA-seq data, and phenotype data (in the form of HPO terms) to find Mendelian disease diganoses caused by splicing mutations in the cohort. This framework can be applied to any dataset (1) that has Mendelian disease patients and (2) has genomic, long-read transcriptomic, and clinical phenotyping data. Because our BDB dataset is enriched in Mendelian disease patients, it is an ideal dataset to apply this framework to.
Please refer to the genetic_disease_diagnosis folder for full details.
We were lucky to have access to a powerful computational cluster to perform these analyses on our dataset. Some of the analyses require a significant level of computational power to complete. Some of our scripts may not be optimized for computational efficiency, and we may improve them in the future.
For questions/concerns/bugs feel free to open an issue in this repository. We will get back to you as soon as possible.
If you use any of our code/methods/analysis in your work, please cite: [INSERT PAPER LINK HERE]
MIT License
Copyright (c) 2025 Children's Hospital of Philadelphia.
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.