RNA-seq analysis

Raw reads from mRNA-seq were first trimmed using fastp with default settings (v0.23.2; (64)). Clean reads were locally aligned to the L. major Friedlin strain genome version 9.0 using bowtie2 tool (65) with ‘very high sensitivity’ parameter and further processed with samtools version 1.7 (66). To ensure proper read placement, alignments with multiple low-quality hits and mapping quality (MQ) scores less than 30 were removed. Strand-specific read coverage was calculated from BAM files obtained from Bowtie2 using customized pysam scripts (https://github.com/pysam-developers/pysam). To compare read coverage for different sites the mean read coverage was calculated and normalized by total number of read number of sequencing library for each category. To compare read coverage of dSSRs, we analyzed the SSRs and the 5-kb flanking regions with DeepTools (3.5.0) using 100 bp bins flanking the SSRs and dividing each SSR into 50 equally sized bins (67). For each sample, FeatureCounts (v2.0.1) (68) was used to count reads for each reference transcript annotation, followed by normalization/variance stabilization using DESeq2 (v1.26.0) (69). Differentially expressed genes (DEGs) were identified using the DESeq2 by comparing WT and PP1 KO samples in triplicate by setting log2 fold change >1 and FDR-adjusted P-value <0.001 (Supplementary Table S1). The average mapping rate for mRNA-seq replicates of WT and PP1 KO was 0.92 (Supplementary Table S2). DMOG RNA-seq data shown here are from previously published work (36).