TFs_predict predicts candidate transcription factors in non-model organisms using transcription-factor annotations from model organisms. The pipeline involves the AnimalTFDB4.0 database and 1-way Best Hits (DIAMOND).

Inputs:

  • Query protein FASTA (all genes) after using extract_isoform pipelines
  • Reference protein FASTA of TFs (from AnimalTFDB4.0)

Workflow:

  1. Download human TF FASTA from AnimalTFDB4.0
  2. Predict TFs of query species fasta using AnimalTFDB4.0 website tool
  3. Find the 1-way Best Hits (DIAMOND) between query species fasta and human TF FASTA (Database): output all query species sequences
  4. Union/Intersection the output of step 2 and 3:
    • Union: includes more potential true TFs, but also more false positives.
    • Intersection: reduces false positives, but may miss more true TFs.

AnimalTFDB 4.0 background

TF prediction on the AnimalTFDB 4.0 website: the database authors’ pipeline

Classification rules described by the AnimalTFDB 4.0 authors

The following background describes the database authors’ method. The analysis here uses their predictions and combines them with a best-hit search.

  • **First, the database authors classified the TFs based on the family-specific domain when a superfamily has several families. For example, for the zf-C2H2 superfamily, the database authors classified proteins with both zf_C2H2 and ZBTB domains into the ZBTB family and proteins containing only zf_C2H2 domain into the zf_C2H2 family. **
  • **The second rule is that if a TF has several different DBDs, it is assigned to the family with the smallest E value of DBD. **
  • The third rule is removing enzyme proteins since they were annotated as enzymes although they contain some DBDs.

AnimalTFDB4.0 Pipelines

Database methods and supplementary materials: AnimalTFDB 4.0.

The database authors built a prediction pipeline using TF families and classification rules

  • **Firstly, the database authors downloaded the Hidden Markov Model (HMM) profiles for DBDs of 58 TF families from the latest Pfam database (version 35.0) (23). ** Each TF has at least one DNA-binding domain (DBD) which is conserved in evolution. Based on their DBDs, TFs could be classified into different families. After reviewing literatures, the database authors collected and curated 72 animal TF families and a group of them named “Others” which includes some orphan TFs. The database authors identified TFs based on the Hidden Markov Model (HMM) profiles of their DBDs. Among the 72 defined families, 59 families had HMM profiles of their DBDs in Pfam database (v31.0) and the database authors downloaded them directly.
  • **Secondly, the database authors reconstructed the remaining 14 TF families based on the DBD sequences from classical species (human, mouse, zebrafish, and fly) with HMMER (v3.1b2) (42). ** For the remaining domains without available Pfam HMM profiles, the database authors rebuilt the HMM profiles using the sequences in representative species (human, mouse, zebrafish and fly). To build the HMM profiles for them, the database authors performed multiple sequence alignment by ClustalW2 for their DBD sequences and used the hmmbuild program in HMMER 3.1b2 package to build HMM profiles.
  • Third, the website primarily runs the HMMER hmmsearch step, searching protein sequences against the DNA-binding-domain HMM profiles. The database authors report family-specific E-value thresholds in their supplementary table. They adjusted some thresholds after updating domain profiles and protein sequences: Fork_head from 1e-4 to 1e-3, and Homeobox to 1e-2. This description concerns the AnimalTFDB authors’ implementation, rather than a new method developed here.
  • Furthermore, the database authors examined the ‘Others’ family and removed some proteins that are not TFs, such as centromere proteins (CENPA, CENPS, CENPT, CENPX).
AnimalTFDB 4.0 authors' diagram of TF and transcriptional-cofactor prediction and annotation
AnimalTFDB 4.0 authors' diagram of TF and transcriptional-cofactor prediction and annotation

The workflow, data and annotation summaries of AnimalTFDB 4.0. The top of the figure illustrates the TFs prediction workflow, and the bottom of the figure describes the TcoFs prediction workflow. The left of the figure shows the data statistics and tools in AnimalTFDB 4.0. The right of the figure shows the multiple types of annotations in AnimalTFDB 4.0, of which variations, gene expression, PTMs and autophagy are newly added.

Example: Bfl

Download human TF FASTA from AnimalTFDB4.0

Predict TFs of query species fasta using AnimalTFDB4.0 website tool

The recorded AnimalTFDB 4.0 submission limit was 10,000 sequences, so split Bfl_proteins_output_final.fa into three FASTA files.

awk 'BEGIN {n=0;} /^>/ {if(n%10000==0){file=sprintf("chunk%d.fa",n);} print >> file; n++; next;} { print >> file; }' < Bfl_proteins_output_final.fa

Run the website prediction three times, copy the results into a spreadsheet, save them as text, and upload them to the analysis server.

~//Bfl_AnimalTFDB4.0_TF_list*
for file in Bfl_AnimalTFDB4.0_TF_list1.txt Bfl_AnimalTFDB4.0_TF_list2.txt Bfl_AnimalTFDB4.0_TF_list3.txt; do
cat $file;
echo;
done > Bfl_AnimalTFDB4.0_TF_list.txt

At least one candidate TF domain was found in 872 + 793 + 358 = 2,023 sequences.

Find the 1-way Best Hits (DIAMOND) between query species fasta and human TF FASTA (Database): output all query species sequences

conda activate //DIAMOND
diamond makedb --in Homo_only_TF_protein_output_final.fasta -d Homo_TF &&

diamond blastp -q Bfl_proteins_output_final.fa -d Homo_TF -o Bfl-vs-HomoTF.tsv --ultra-sensitive --matrix BLOSUM90 --iterate

For each Bfl sequence I only select the best hit base on bit score

awk '{print $1}' Bfl-vs-HomoTF.tsv | sort -u | while read -r gene_id; do

best_hit=$(grep "\\b${gene_id}" Bfl-vs-HomoTF.tsv | awk -F "\t" '$11<=1e-4' | awk -F "\t" '$3>=30' | sort -k 12n | tail -1 | awk '{print $1, $2}' OFS="\t")

echo "$best_hit"

done > Bfl_BH.txt


grep . Bfl_BH.txt > Bfl_BH_rmBlank.txt
awk 'END {print NR}' Bfl_BH_rmBlank.txt

→ 1863 genes

Union of the one-way best hits and Bfl_AnimalTFDB4.0_TF_list.txt

venn union <(awk '{print $1}' Bfl_AnimalTFDB4.0_TF_list.txt) <(awk '{print $1}' Bfl_BH_rmBlank.txt) > Union_Bfl.AnimalTFDB4.0_Bfl.BH.txt

awk 'END {print NR}' Union_Bfl.AnimalTFDB4.0_Bfl.BH.txt

→ 2509 genes Compare with the available Bfl TF list of 1,928 genes, treating that list as a reference for this comparison.

venn intersection <(grep ">" ~/Shared/R145_SchneiderSQ/TO-GCN/TF-OrthoFinder/processed_protein_fasta/putative_Bfl_TF_output_1.9.fa)  <(awk '{{ $0=">"$0 } print}' Union_Bfl.AnimalTFDB4.0_Bfl.BH.txt) | wc -l
  • 1,924 genes overlap: 1,924 / 2,509 = 76.68% of the predicted union is in the reference list.
  • 585 genes occur only in Union_GRCh38.AnimalTFDB4.0_GRCh38.HomoTF.RBH.txt.
  • Four genes occur only in the reference Bfl TF list.

Intersection of the one-way best hits and Bfl_AnimalTFDB4.0_TF_list.txt

venn intersection <(awk '{print $1}' Bfl_AnimalTFDB4.0_TF_list.txt) <(awk '{print $1}' Bfl_BH_rmBlank.txt ) > Intersection_Bfl.AnimalTFDB4.0_Bfl.BH.txt

awk 'END {print NR}' Intersection_Bfl.AnimalTFDB4.0_Bfl.BH.txt

→ 1377 genes Compare with the available Bfl TF list of 1,928 genes, treating that list as a reference for this comparison.

venn intersection <(grep ">" ~/Shared/R145_SchneiderSQ/TO-GCN/TF-OrthoFinder/processed_protein_fasta/putative_Bfl_TF_output_1.9.fa)  <(awk '{{ $0=">"$0 } print}' Intersection_Bfl.AnimalTFDB4.0_Bfl.BH.txt) | wc -l
  • 1,365 genes overlap: 1,365 / 1,377 = 99.13% of the predicted intersection is in the reference list.
  • 12 genes occur only in Intersection_GRCh38.AnimalTFDB4.0_GRCh38.HomoTF.RBH.txt.
  • 563 genes occur only in the reference Bfl TF list.
← Back to projects