https://doi.org/10.5061/dryad.0k6djhbcm

This repository contains the data for the paper:

"Illuminating the mystery of thylacine extinction: a role for relaxed selection and gene loss"

DOI: 10.1098/rspb.2025.1339

Folder Structure

• Main folder: Thylacine_Gene_Loss.zip
• Chromosome_wise_Chains_and_TOGA_results: Contains the output files of make_lastz_chains and TOGA.
• Purpose: Contains input, intermediate, and final output files from genome alignment and ortholog projection using LASTZ and TOGA.
• Key Tools:
• make_lastz_chains (https://github.com/hillerlab/make_lastz_chains): pairwise genome alignment
• TOGA (https://github.com/hillerlab/TOGA): identifies orthologous genes and infers gene status
• Data: 
• Chain files: pairwise genomic alignment between Tasmanian devil and Tasmanian tiger.
• TOGA outputs per chromosome: predicted gene ortholog/gene status, inactivating mutations
• Figure1: Contains the gene loss confirmation genes losses in thylacine and other marsupial species, IGV-reports, and gene loss history. (Figure 1, figure S2, figure S5, figure S7, figure S8 and figure S9)
• Purpose: Stores all data used in Figure 1 and relevant supplementary figures.
• Content:
• IGV screenshots validating gene loss visually
• Gene loss summary tables
• Gene loss history across lineages
• Figure2: Contains the figure 2 figures and the inputs required to get figures like phylogenetic trees from TimeTree and iTOLs. (Figure 2)
• Content:
• Phylogenetic trees generated using the TimeTree website (www.timetree.org).
• Tree annotations and visualisations with iTOL (https://itol.embl.de/).
• Figure3:
• (a) Contains the input and output files of GC-content variation in thylacine. (Figure 3a)
• Analyses genome-wide GC content of genes, especially in lost genes and gene status inferred using TOGA.
• Input: gene coordinates and sequences
• Output: scatter plot, box plots, violin plots (GC% V/s GC stretches per gene), and codes associated with these plots.
• (b) Selection pressure analysis of "clearly lost" genes (found using TOGA), tested using RELAX. (Figure 3b)
• Uses RELAX to test for relaxation or intensification selection
• Input: codon alignments (built using TOGA) and gene trees (obtained from the TimeTree website).
• Output: Test statistics (k, p-values), model fits, output files.
• (c) Expression of lost genes in Devil and Dunnart transcriptome. (Figure 3c)
• Compares expression of lost genes in Devil and Dunnart using RNA-seq data.
• Input: Gene read counts (from kallisto), normalised expression.
• Output: Heatmaps showing tissue-specific expression.
• Gene_Tree_SAMD9-9L: Contains the alignment and tree files used to obtain gene trees. (figure S3)
• Hybpiper_patchwork: Contains the input files for Hybpiper and patchwork. The subset of fastq reads gives hits (BLASTn) to lost genes (as query sequences).
• Hybpiper_SAMD9-9L: In-frame stop codon in SAMD9L reassessed using short-read data of thylacine, using SAMD9 and SAMD9L coding sequences of other 20 marsupial species.
• Purpose: Re-evaluate premature stop codons in SAMD9L using thylacine short-read data and 20 marsupial orthologs
• Content:
• Alignment of CDS from other marsupials.
• Read-based confirmation of frameshifts or nonsense mutations.
• IGV_reports-BWA+BLASTn: Contains the IGV-reports for clearly lost genes in chromosome-wise folders (Figure 1)
• Purpose: Visualisation of "clearly lost" genes using Integrated Genomics Viewer (IGV)
• Input: Short-read BAM files (built using blastn 17 output format or mapped with bwa), genome references of Tasmanian devils.
• Figures: Visual validation of frame-disrupting changes shown in Figure 1.
• Loss_confirmation: Contains the output files of short-read BLASTn, bam-read count, and gene-loss events found by TOGA to generate a summary table and confirm gene loss. (Figure 1, table S1)
• Purpose: Independent confirmation of gene losses detected by TOGA
• Components:
• bam-readcount results: coverage and mutation confirmation
• BLASTn alignments: confirm absence/presence of genes
• Summary table: includes gene ID, status, evidence, and references
• Main_Figures: contains the biorender permission for publication.
• Myrmecobius_fasciatus_RNA_seq: Contains an IGV screenshot for the VWA7 gene, subseted bam files, and BED files. (figure S10)
• PAML_omega_SAMD9-9L: Evolutionary rate dN/dS ratio (omega) across the SAMD9 and SAMD9L genes. Contains the input and R script. (figure S4)
• Purpose: Molecular evolution analysis of SAMD9/9L using PAML
• Output:
• ω (dN/dS) values
• Selection test comparisons
• Scripts: R code used to parse and plot outputs
• Selection_analysis: Selection pressure analysis of lost gene across 21 marsupial species, tested using RELAX, codeML, MEME, FEL, aBSREL, BUSTED, and gBGC. (Figure 3 and table S5)
• Comprehensive selection pressure tests across 21 marsupial species
• Tools Used:
• RELAX, codeML, MEME, FEL, aBSREL, BUSTED, gBGC
• Input: gene trees, codon alignments
• Output: statistical test results, visualisation files
• Thylacine_transcriptome-miRNA: Contains BLASTn search output for lost genes and control genes in the miRNA database of thylacine.
• Purpose: Evaluates whether lost genes are expressed in the transcriptome or not
• Method: BLASTn against thylacine miRNA database
• TOGA_output_sorted: Contains the sorted output files of TOGA, such as inact_mut_data.txt and loss_summ_data.tsv, codon.fasta, nucleotide.fasta.tgz, prot.fasta, query_annotation.bed, and query_gene_spans.bed. Also, a mutation plot was generated using TOGA. (Figure 1)
• Cleaned and sorted outputs from TOGA
• Key files:
• inact_mut_data.txt: details on inactivating mutations
• loss_summ_data.tsv: summary of lost/retained genes (I, M, PI, UL, L, PM)
• codon.fasta, prot.fasta: sequences used for alignment of reference and query species.
• query_annotation.bed: gene annotation intervals
• Includes mutation plots obtained by running TOGA.
• Synteny: Contains the NCBI synteny images for SAMD9L, HSD17B13, CUZD1 and VWA7. (figure S5 and tables S2)
• Screenshots from NCBI Genome Data Viewer (GDV) for gene SAMD9L, CUZD1, HSD17B13, and VWA7.
• Used to identify 1-to-1 gene orthologs.
• SAMD9_gene_loss_in_carnivore: Contains the script for gene loss events identification by TOGA. Assembly verification output files of klumpy and IGV reports of events by short read for dingo and tiger. (figure S6)
• Validates gene loss events for SAMD9 in carnivores
• Contains:
• TOGA gene-loss identification script
• Assembly validation (klumpy)
• IGV-based confirmations
• R1_RSPB: Includes input, output, and code for paralog and GC-content comparison, revalidation of gene loss events, RNA-seq data from Macropus eugenii, and phylogenetic correlation analysis (Figure 2b)
• Supporting analyses for:
• Paralog vs. GC-content trends
• RNA-seq revalidation from Macropus eugenii
• Phylogenetic correlation of gene loss and carnivorous diet.

• Details of tools, version, and purpose

Prerequisites:

• TOGA (1.1.7)
• make_lastz_chains (https://github.com/hillerlab/make_lastz_chains.git)
• PAML (4.9f)
• BLASTN (2.13.0)
• phastBias(1.6)
• mapnh(1.3.0)
• HYPHY 2.5.48(MP) and 2.5.62(MP)
• BWA(0.7.17-r1188)
• Samtools (1.3)
• bedtools (v2.27.1)
• ea-utils (https://github.com/ExpressionAnalysis/ea-utils.git)
• seqtk (1.2-r94)
• Guidance (v2.01)
• PRANK (v.170427)
• IGV (2.8.13)
• igv-reports (1.12.0)
• STAR (2.7.0d)
• MegaX
• R (4.1.0)
• bam-readcount (1.0.1)
• HybPiper (2.3.1)
• Patchwork (0.5.2)
• IQ-TREE (2.3.6)
• klumpy (1.0.11)
• kallisto (0.51.0)
• Bowtie2 (2.5.4)

R packages:

• ape (5.5)
• phytools (0.7-90)
• ggplot2 (3.3.5)
• ggrepel (0.9.1)
• cowplot (1.1.1)
• dplyr (1.0.7)
• ggplotify (0.1.0)
• grid (4.1.1)
• gridExtra (2.3)
• reshape2 (1.4.4)