CDS Database and Custom Scripts for a genomics paper: Pangolin genomes offer key insights and resources for the conservation of the world's most trafficked mammal (submitted)
Creators: name: Heighton, Sean Paul
name: Allio, Rémi
name: Murienne, Jérôme
name: Salmona, Jordi
name: Meng, Hao
name: Scornavacca, Celine
name: Bastos, Armanda D.S
name: Njiokou, Flobert
name: Pietersen, Darren William
name: Tilak, Marie-Ka
name: Luo, Shu-Jin
name: Delsuc, Frédéric
name: Gaubert, Philippe
Publisher: name: Zenodo
Year: 2022
Subjects: Pangolin
Pholidota
genomics
phylogeny
full genomes
genomic diversity
conservation
CDS
PSMC
introgression
divergence estimates
forensic markers
wildlife trade

description: As part of our publication on the conservation genomics of pangolins (submitted), we provide a database of CDSs ranked by diversity amongst the eight pangolin species and set of custom scripts created and implemented during the processing of the genomic data. More information on when and how these scripts were implemented can be found in the methods section or the pipeline outlines in Figure S1. <strong>Database S1:</strong> An excel file with a list of 11 741 Coding Sequences (CDS) of Genes ranked from most diverse to least diverse amongst the eight pangolin species (Pholidota) based on mean pairwise identities (sheet "CDSs"). <em>PhyKIT</em> and <em>DnaSP v6</em> outputs were merged based on CDS Gene ID (Ensembl) with a custom script (Custom script 4) to obtain a range of diversity statistics per CDS marker. Outlying genes with low levels of mean pairwise identity (or high levels of diversity between species) were removed (Figure S2) and stored in a separate sheet within the Database called "Removed CDSs". Please refer to the methodology section of the manuscript for more information on the methods used to obtain this database. Column headings are based on <em>PhyKIT</em> (https://github.com/JLSteenwyk/PhyKIT) and <em>DnaSP v6 </em>(http://www.ub.edu/dnasp/) outputs, more information on these statistics and how they are estimated can be found in their respective manuals. <strong>Custom script 1:</strong> Removes paired reads that matched a different chromosome/scaffold (CDS marker in our case) by keeping only the records where column 7 in the sam/bam file contained an '=' (mate matches mapped on the same CDS) with an <em>awk</em> command. This required extracting the header with a <em>grep</em> command before and then attaching it again with a <em>cat</em> command at the end (Figure S1a). <strong>Custom script 2:</strong> Obtain the average depth per CDS marker (gene/scaffold) and the percent of CDS marker covered at a certain minimum or maximum depth of your choosing. We chose to identify the percent of each CDSs covered by 10x depth or more. This will provide an output of CDS id's and their respective values of depth and percent coverage at a certain depth. It requires the use of <em>samtools depth</em> before running this script in order to obtain the depth at each nucleotide position of each CDS (Figure S1a) <strong>Custom script 3:</strong> This provides a list of CDSs (genes/scaffolds) that have too high/ too low average mapping depth of your choice (indicative of repetitive regions, bad reference sequence, paralogs etc) and too little percent of the CDS covered at a depth of your choice. It then combines the two lists, removes duplicates by CDS id and then uses <em>seqkit</em> to remove these CDSs from the species consensus fasta file comntaining all the CDSs for a specific species. The choice of data relies on the custom script 2 outputs and your intuition on what should be the cutoffs (we used a maximum average depth cut-off 2x the mean depth for all CDSs for each species for example; Figure S1a). <strong>Custom script 4:</strong> This provides a list of CDSs (genes/scaffolds) with all relevant information to indicate CDS diversity and selection pressures extracted from <em>phykit</em> and <em>DNAsp v6</em> outputs (Figure S1a). The list can then be sorted by whichever diversity or selection index you prefer (our list can be viewed in the Database S1 containing). <strong>Custom script 5:</strong> Obtain the average depth of the genome mapping and the percent of genome covered at a certain minimum or maximum depth of your choosing. This is the same as custom script 2 but is for the entire genome rather than for individuals CDS markers (Figure S1b). It also requires the use of <em>samtools depth</em> before running this script in order to obtain the depth at each nucleotide position of the mapped genome. Additionally it is in the form of a loop script for multiple samples to be run in parallel. 5a calls the loop using the id.txt of all the samples you're interested in running in 5b. If the loop is not required, you can run 5b buy changing 'ind=$1' to 'ind=YOUR_SAMPLE'.