Motif activity calculation in 3 steps. -------------------------------------- These scripts require a recent version of Python (later than 3.0), as well as a recent version of Numerical Python (numpy). Input files ----------- promoters.bed: BED file containing the position of the transcription start site of each gene. For motif activity calculations on enhancers instead of promoters, the BED file also must contain the midpoint position of the enhancer in the thickStart and thickEnd columns (columns 7 and 8), as generated by Robin Andersson's enhancer prediction scripts (see Andersson et al.: "An atlas of active enhancers across human cell types and tissues", Nature 507 (7493): 455-461 (2014)). expression.txt: Tab-delimited file containing the expression data. All lines in this file must have the same number of columns. The first column of each line must specify the promoter name, which must be consistent with the gene names used in the promoters.bed file. The text in the first column of the first line is arbitrary; the remaining columns on the first line must contain the sample names. tfbs.bed.bz2: File containing the precalculated TFBSs for your genome assembly; the file must be bzipped BED file. Precalculated TFBSs for human genome hg19 and hg38, mouse genome mm9 and mm10, rat genome rn6, dog genome canFam3, and chicken genome galGal5 and galGal6 are available as part of the Supplementary Materials of Alam et al., "Comparative transcriptomics of primary cells in vertebrates", Genome Research 30(7): 951-961 (2020). Step 1: Calculate the binding profile of TFs with respect to the TSS -------------------------------------------------------------------- Use the script make_profile.py to calculate for each TF the binding profile with respect to the transcription start site: python make_profile.py --output=profile.txt --input=promoters.bed --tfbs=.sites.bed.bz2 where promoters.bed is the BED file containing the position of the transcription start site of each gene, and .sites.bed.bz2 is a BED file with the genomic locations of the predicted TFBSs. The script creates a tab-delimited file containing the binding profiles for all transcription factors with predicted binding sites in .sites.bed.bz2. The full list of options is: python make_profile.py -t [-u ] [-d ] -i -o [-c cutoff] [-p] [-s] Options: --tfbs= or -t Path to the file containing the precalculated TFBSs for your genome assembly; the file should be bzipped BED file. --upstream= or -u Distance upstream of the TSS to include in the profile. Defaults to 300. --downstream= or -d Distance downstream of the TSS to include in the profile. Defaults to 100. --input= or -i Name of the input file containing the transcription start site positions in BED file format. --output= or -o Name of the output file in which the calculated profiles should be stored. --cutoff= or -c Required minimum posterior probability of a predicted transcription factor binding site. Predicted transcription factor binding sites with a score (as defined in the BED file) lower than this cutoff will be ignored. Default value is 0. --symmetric or -s Calculate a symmetric profile. This requires upstream and downstream to be equal to each other. --plot or -p Create a figure for each TF showing its binding density. For enhancers, a symmetric profile from -500..+500 base pairs can be calculated using python make_profile.py --output=profile.txt --input=enhancers.bed --tfbs=.sites.bed.bz2 --upstream=500 --downstream=500 --symmetric Step 2: Associate TFBSs with promoters -------------------------------------- Use the script associate_tfbs.py to associate predicted TFBSs to CAGE promoters: python associate_tfbs.py --output=tfbs.txt --input=promoters.bed --tfbs=.sites.bed.bz2 --profile=profile.txt where promoters.bed is the BED file containing the position of the transcription start site of each gene, .sites.bed.bz2 is a BED file with the genomic locations of the predicted TFBSs, and profile.txt is the TFBS profile calculated in Step 1. The window around the transcription start site to search for transcription factor binding sites is taken from the profile.txt file. The full list of options is: python associate_tfbs.py -t -p [-c cutoff] -i -o Options: --tfbs= or -t Path to the file containing the precalculated TFBSs for your genome assembly; the file should be bzipped BED file. --profile= or -p Name of the file in which the precalculated binding profiles with respect to the CAGE-defined TSS are stored. --cutoff= or -c Required minimum score of a predicted transcription factor binding site. Predicted transcription factor binding sites with a score lower than this cutoff will be ignored. Default value is 0. --input= or -i Name of the input file containing the transcription start site positions in BED file format. --output= or -o Name of the output file in which the associated TFBSs will be stored. Step 3: Calculate the motif activities -------------------------------------- We are now ready to calculate the motif activities: python calculate_motif_activity.py --output=activities.txt --input=expression.txt --tfbs=tfbs.txt where expression.txt is the tab-delimited file containing the expression data, and tfbs.txt is the name of the file created in Step 2. This script will generate the file activities.txt containing the motif activities in each condition. The full list of options is: python calculate_motif_activity.py -i -o -t [-n ] [-q] Options: --t= or -t Name of the file containing the predicted TFBSs for each promoter. --input= or -i Name of the input file containing the expression level for each promoter in tags-per-million [tpm]. --output= or -o Name of the output file in which the calculated profiles should be stored. --n= or -n Minimum number of predicted TFBSs for each motif. Motifs with fewer predicted TFBSs than this number are ignored. Default value is 150.