Introduction
============

One of the most ancient and abundant classes of DNA binding domains (DBDs) is the C2H2 zinc finger class ([@evv141-B97]; [@evv141-B92]; [@evv141-B28]).  The C2H2 zinc finger domain has two cysteine and two histidine residues that coordinate a zinc ion, and typically consists of the amino acid sequence C-X(2-4)-C-X(12)-H-X(3-5)-H ([@evv141-B15]; [@evv141-B73]).  C2H2 zinc finger motifs are found in many transcription factors and based on their arrangement and number can be subdivided into different families ([@evv141-B51]).  The Krüppel-like factor and specificity protein (KLF/SP) transcription factor gene family is characterized by a highly conserved triple-C2H2 DBD located toward the C-terminus composed of three tandem zinc fingers that are evenly spaced by conserved linker regions ([@evv141-B51]) and share similarity with the *Drosophila Krüppel* gene ([@evv141-B96]).  This C2H2 zinc finger DBD (KLF-DBD) binds to Guanine-Cytosine-rich regions and CACC elements (GT boxes) ([@evv141-B56]).  The more N-terminal regions of KLF/SP transcription factors are typically highly variable and consist of different combinations of transactivation/repression domains.  Historically, mammalian KLFs have been divided into 3 groups based on shared domain architecture: The KLF1, 2, 4, 5, 6, and 7 groups; the KLF3, 8, and 12 groups; and the KLF9, 10, 11, 13, 14, 16 groups ([@evv141-B70]), whereas SPs, which differ from KLFs by the presence of the Buttonhead (Btd) box domain just 5′ of the KLF-DBD, are typically divided into 2 groups: SP1--4 and SP5--9 ([@evv141-B123]). 

KLF/SP genes within each domain architecture group share similar functions based on the retention of explicit transactivation motif complements.  A range of studies present a complex picture in which KLF/SP genes can be singly or combinatorially involved in temporally and spatially disparate cellular and developmental processes.  For example, fly embryos mutant for *luna*, the *Drosophila* KLF6/7 ortholog, die early during development due to mitotic defects ([@evv141-B27]; [@evv141-B130]), whereas *cabut*, a KLF9/13 ortholog in *Drosophila*, is required for proper dorsal closure during gastrulation ([@evv141-B75]).  However, *cabut* also plays a role later in fly organ development by coordinating signaling for proper wing disc patterning ([@evv141-B94]).  Among the vertebrates, KLF genes are often associated with balancing stem cell proliferation and differentiation, as well as regulating metabolic homeostasis.  The most notable member is KLF4, one of the four pioneer transcription factors required to induce pluripotency in human and mouse fibroblasts ([@evv141-B125]; [@evv141-B116], [@evv141-B117]) and a component of a core circuit of genes that maintain self-renewal in mammalian embryonic stem cells along with KLF2 and KLF5 ([@evv141-B52]).  However, in gut epithelia, KLF4 regulates terminally differentiated cells while KLF5 is expressed in the proliferating crypt cells ([@evv141-B69]).  In mammals, KLF2 together with KLF1 and KLF13 also regulate erythrocyte maturation and differentiation as well as globin gene activity ([@evv141-B72]; [@evv141-B10]; [@evv141-B41]).  KLF2 in zebrafish contributes to the differentiation of ectoderm derived tissues ([@evv141-B61]).  In mammals, including humans, KLF11 and KLF14 play an important role in the regulation of genes associated with diabetes and metabolic syndrome phenotypes, respectively ([@evv141-B110]; [@evv141-B65]).  Similarly, complex intersections with both development and metabolism exist for members of the SP subfamily.  For example, in mammals, SP1, SP3, and SP7 regulate osteoblast mineralization and differentiation ([@evv141-B76]; [@evv141-B124]).  SP1 is also an important regulator of metabolic genes involved in the glycolytic pathway, fatty acid synthesis, and ribosome biogenesis ([@evv141-B7]; [@evv141-B78]).  Overall, members of the KLF/SP gene family are known to function in a wide variety of biological processes ([@evv141-B12]; [@evv141-B137]; [@evv141-B82]; [@evv141-B131]; [@evv141-B70]; [@evv141-B138]; [@evv141-B126]). 

In contrast to the extensive studies highlighting the importance of the KLF/SP genes to core cellular processes, comparatively few studies have investigated the evolutionary relationship of KLF/SP genes in lineages outside of mammals ([@evv141-B60]; [@evv141-B66]; [@evv141-B105]; [@evv141-B17]; [@evv141-B71]; [@evv141-B100]; [@evv141-B103]).  A KLF gene was recently identified in the choanoflagellate *Monosiga brevicollis* genome; however, that study's conclusions were restricted to examining porcine KLF paralogy ([@evv141-B17]).  A more recent study, focused on the phylogenetic distribution of C2H2 zinc finger families in eukaryotes, also showed that KLFs were present in *Monosiga* but absent in the fungal taxa surveyed ([@evv141-B104]).  No study to date has examined the phylogenetic context of the different transactivation/repression domains associated with the KLF/SP gene family.  Pinpointing the origin and evolutionary history of this gene family and associated domains can help determine possible relationships of the KLF/SP repertoire expansion to key innovations in the evolution of metazoan cellular diversity.  Hypotheses of metazoan gene evolution are greatly aided by sampling a wide range of taxa that include nonmetazoan representatives.  Here, we infer the evolutionary history of the KLF/SP gene family and their associated transactivation/repression domains across a wide range of eukaryotes.  Species were chosen to represent well established groups across the Eukarya including the following: Bikonta, Amorphea (Amoebozoa + Apusozoa + Opisthokonta), Opisthokonta (Holomycota + Holozoa), and Holozoa (Ichthyosporea + Filasterea + Choanoflagellata + Metazoa) ([fig.  1](#evv141-F1){ref-type="fig"}) ([@evv141-B3]; [@evv141-B30]; [@evv141-B80]).  A number of transactivation/repression domains found in KLF/SPs are also present in the genomes of ancient unicellular lineages that lack KLF/SPs, lending support for domain shuffling playing a major role in the acquisition of transactivation/repression domains during the expansion and diversification of the KLF/SP gene family.  We show that domain connectivity and resulting unique domain architectures among these transcription factors have become increasingly complex in metazoans.  Thus a pattern of gene duplication along with domain shuffling and the rare emergence of de novo domains have collectively played a vital role in the evolution of the KLF/SP gene family.  F[ig]{.smallcaps}.  1.---Distribution of C2H2 zinc finger proteins, KLF-DBD containing proteins, and KLF/SP proteins in representative Eukarya taxa.  Rows indicate representative genomes searched.  Columns indicate the total number of protein sequences that contain at least one C2H2 zinc finger using the Pfam PF00096 HMM model, the total number of protein sequences that contain the archetypical KLF-DBD, the total number of bona fide KLF sequences recovered, and the total number of SP sequences recovered.  Phylogeny is based on [@evv141-B3], [@evv141-B30], [@evv141-B31], [@evv141-B99], and [@evv141-B102]. 

Materials and Methods
=====================

KLF Identification Pipeline
---------------------------

To broadly identify C2H2 zinc finger and KLF/SP proteins associated with a diverse range of Eukaryote lineages, publicly available sequenced genomes listed in [table 1](#evv141-T1){ref-type="table"} were comprehensively searched including 26 metazoans, 4 unicellular holozoans, 7 holomycotans, 1 apusozoan, 4 amoebozoans, and 6 bikonts.  We used the HMMER 3.0 program ([@evv141-B32]), to identify proteins that contained C2H2 zinc fingers.  The Hidden Markov Model (*hmm*) of the C2H2 zinc finger domain PF00096 ([@evv141-B87]) was downloaded from the Pfam database.  The *hmmsearch* command was then used to search protein models of the representative 48 Eukaryote species using the PF00096 *hmm* as a query.  HMMER identified all protein sequences (using default settings) that contained at least one C2H2 zinc finger corresponding to the *hmm.*  This output was used for subsequent analyses.  Raw outputs for this and subsequent steps of the pipeline can be found in [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online Table 1Species used in this Study with Reference to Genome or Transcriptome DatabaseSpeciesGenome/TranscriptomeReferenceAmorphea*Homo sapiens*<ftp://ftp.ncbi.nlm.nih.gov/genomes/H_sapiens/>[@evv141-B63]*Mus musculus*[www.ensembl.org/Mus_musculus/](http://www.ensembl.org/Mus_musculus/)[@evv141-B18]*Xenopus tropicalis*<http://useast.ensembl.org/Xenopus_tropicalis/Info/Index>[@evv141-B46]*Latimeria chalumnae*<http://useast.ensembl.org/Latimeria_chalumnae/Info/Index>[@evv141-B4]*Danio rerio*[www.ensembl.org/Danio_rerio/Info/Index](http://www.ensembl.org/Danio_rerio/Info/Index)[@evv141-B48]*Takifugu rubripes*[www.ensembl.org/Takifugu_rubripes/Info/Index](http://www.ensembl.org/Takifugu_rubripes/Info/Index)[@evv141-B5]*Gasterosteus aculeatus*<http://useast.ensembl.org/Gasterosteus_aculeatus/Info/Index>[@evv141-B53]*Callorhinchus milii*<http://esharkgenome.imcb.a-star.edu.sg/download/>[@evv141-B128]*Leucoraja erinacea*<http://skatebase.org/downloads>*Scyliorhinus canicula*<http://skatebase.org/downloads>*Petromyzon marinus*<http://useast.ensembl.org/Petromyzon_marinus/Info/Index>[@evv141-B113]*Ciona intestinalis*<http://useast.ensembl.org/Ciona_intestinalis/Info/Index>[@evv141-B29]*Branchiostoma floridae*<http://genome.jgi-psf.org/Brafl1/Brafl1.home.html>[@evv141-B88]*Strongylocentrotus purpuratus*<http://metazoa.ensembl.org/Strongylocentrotus_purpuratus/Info/Index>[@evv141-B114]*Daphnia pulex*<http://genome.jgi-psf.org/Dappu1/Dappu1.download.ftp.html>[@evv141-B20]*Drosophila melanogaster*<http://flybase.org>[@evv141-B2]*Tribolium castaneum*<http://beetlebase.org>[@evv141-B93]*Capitella teleta*<http://genome.jgi.doe.gov/Capca1/Capca1.download.ftp.html>[@evv141-B109]*Lottia gigantea*<http://genome.jgi-psf.org/Lotgi1/Lotgi1.download.ftp.html>[@evv141-B109]*Nematostella vectensis*<http://genome.jgi.doe.gov/Nemve1/Nemve1.home.html>[@evv141-B89]*Acropora digitifera*<http://www.compagen.org/datasets.html>[@evv141-B106]*Trichoplax adhaerens*<http://genome.jgi.doe.gov/Triad1/Triad1.home.html>[@evv141-B118]*Amphimiedon queenslandica*<http://spongezome.metazome.net/cgi-bin/gbrowse/amphimedon/>[@evv141-B119]*Oscarella carmela*<http://www.compagen.org/datasets.html>*Mnemiopsis leidyi*<http://research.nhgri.nih.gov/mnemiopsis/>[@evv141-B99]*Pleurobrachia bachei*<http://neurobase.rc.ufl.edu/pleurobrachia>[@evv141-B74]Metazoa*Monosiga brevicollis*<http://genome.jgi.doe.gov/Monbr1/Monbr1.home.html>[@evv141-B59]*Salpingoeca rosetta*<http://www.broadinstitute.org/annotation/genome/multicellularity_project/MultiDownloads.html>[@evv141-B34]*Capsaspora owczarzaki*<http://www.broadinstitute.org/annotation/genome/multicellularity_project/MultiDownloads.html>[@evv141-B122]*Sphaeroforma arctica*<http://www.broadinstitute.org/annotation/genome/multicellularity_project/MultiDownloads.html>Holomycota*Ustilago maydis*<http://www.broadinstitute.org/annotation/genome/ustilago_maydis>[@evv141-B57]*Saccharomyces cerevisiae*<http://www.yeastgenome.org/>[@evv141-B40]*Mortierella verticillata*<http://www.broadinstitute.org/annotation/genome/multicellularity_project/MultiDownloads.html>*Allomyces macrogynus*<http://www.broadinstitute.org/annotation/genome/multicellularity_project/MultiDownloads.html>*Spizellomyces punctatus*<http://www.broadinstitute.org/annotation/genome/multicellularity_project/MultiDownloads.html>*Encephalitozoon cuniculi*<http://genome.jgi.doe.gov/pages/dynamicOrganismDownload.jsf?organism=Enccu1>[@evv141-B58]*Fonticula alba*<http://www.broadinstitute.org/annotation/genome/multicellularity_project/MultiDownloads.html>Opisthokonta*Thecamonas trahens*<http://www.broadinstitute.org/annotation/genome/multicellularity_project/MultiDownloads.html>*Entamoeba histolytica*<http://amoebadb.org/common/downloads/>[@evv141-B64]*Dictyostelium purpureum*<http://genome.jgi-psf.org/Dicpu1/Dicpu1.home.html>[@evv141-B121]*Polysphondylium pallidum*<http://genomes.dictybase.org/pallidum/current>[@evv141-B45]Amorphea*Acanthamoeba castellani*<http://amoebadb.org/common/downloads/>[@evv141-B19]Bikonta*Naegleria gruberi*<http://genome.jgi-psf.org/Naegr1/Naegr1.download.ftp.html>[@evv141-B37]*Cyanidioschyzon merolae*<http://merolae.biol.s.u-tokyo.ac.jp/download/>[@evv141-B67]*Arabidopsis thaliana*<http://www.arabidopsis.org/>[@evv141-B6]*Paramecium tetraurelia*<http://paramecium.cgm.cnrs-gif.fr/download/fasta/>[@evv141-B9]*Bigelowiella natans*<http://genome.jgi.doe.gov/Bigna1/Bigna1.download.ftp.html>[@evv141-B25]*Phytophtora infestans*<http://protists.ensembl.org/Phytophthora_infestans/Info/Index>[@evv141-B42]

We then used a perl script modified from [@evv141-B135] to search the protein sequences of the HMMER output for the 81 amino acid triple-C2H2 zinc finger DBD conserved in KLF/SPs (KLF-DBD).  The amino acid sequence is as follows: C-X~4~-C-X~12~-H-X~3~-H-X~7~-C-X~4~-C-X~12~-H-X~3~-H-X~7~-C-X~2~-C-X~12~-H-X~3~-H, where X can be any amino acid.  Sequences meeting the following criteria were initially considered putative KLF/SPs: 1) Presence of the KLF-DBD and 2) presence of only three zinc fingers (i.e., no additional zinc fingers other than the KLF-DBD).  Sequences fitting these initial criteria were only found within the Filozoa ([supplementary fig.  S1*A*](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online).  We additionally found 13 nonfilozoan amorphean sequences possessing a KLF-DBD; however, in all cases these sequences also contained additional zinc fingers ([Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online).  To determine the relationship between the nonfilozoan amorphean sequences and the putative filozoan KLF/SPs, all sequences were aligned and included in our phylogenetic analyses. 

Where possible, partial KLF/SP gene models were manually extended by obtaining genomic scaffold regions spanning the location of incompletely annotated gene models.  Briefly, the CLCBio software package was used to map partial KLF/SP nucleotide sequences to genomic scaffolds and identify associated open reading frames.  Extended exonic nucleotide sequences were then translated to obtain enhanced annotation of KLF/SP amino acid sequences.  Protein schematics were created using the PROSITE MyDomians Image Creator (<http://prosite.expasy.org/cgi-bin/prosite/mydomains/>, last accessed April 30, 2014).  Scripts for the KLF/SP identification pipeline are publicly available ([Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online). 

Transactivation/Repression Domain Identification and Characterization
---------------------------------------------------------------------

Custom perl scripts were used to search for KLF/SP-associated transactivation/repression domains.  The motifs used for the perl scripts followed these amino acid sequences: Btd box, C-X-C-P-X-C ([@evv141-B132]); SP box, S/T-P-L-X-ϕ-L-X-X-X-C-X-R/K-ϕ ([@evv141-B44]); SID, D-X~1-4~-X-A-ϕ-X-X-L-M/V/L/A-X-F/M/L/I ([@evv141-B136]); PVDLS, P-V-D-L-S/T ([@evv141-B24]); R2, S-V-I-R-H-T-X-D/E ([@evv141-B22]); and R3, ϕ-X-X-G-X-ϕ-ϕ-ϕ-ϕ-ϕ-P/S-Q/P ([@evv141-B22]), where X can be any amino acid and ϕ represents one of V, I, L, M, F, W, G, A, or P. To confirm the accuracy of the automated identification of transactivation/repression domains, we first searched the completely annotated human KLF/SP amino acid sequences.  In all cases, our searches identified the human KLF/SP sequences that had previously been characterized and annotated as containing the specified domains.  For example, the search with the Btd box model identified all 9 human SP sequences.  Thus, our models are both specific and sensitive enough to capture the intended domains.  To identify the nine-amino-acid transactivation domain (9aaTAD), we utilized a prediction tool from <http://www.med.muni.cz/9aaTAD/index.php> (last accessed April 30, 2014), using the "moderately stringent pattern" option ([@evv141-B85]).  Compositionally biased, low complexity regions (LCRs) associated with KLF/SPs were identified using four independent assessments: The CAST algorithm, <http://athina.biol.uoa.gr/cgi-bin/CAST/cast_cgi> ([@evv141-B86]); ScanProsite, <http://prosite.expasy.org/scanprosite/> (last accessed April 30, 2014) ([@evv141-B26]); the SEG algorithm, <http://mendel.imp.ac.at/METHODS/seg.server.html> (last accessed April 30, 2014) ([@evv141-B133]); and manual curation based on criteria from [@evv141-B108].  A putative LCR was assigned if it was identified in at least two of the four methods ([supplementary table S2](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online). 

Additionally, we searched the entire genomes from our representative 48 eukaryotes for the following transactivation/repression domains in proteins excluding identified KLF/SP genes: Btd box, SP box, SID, PVDLS, R2, and R3 domain motifs.  Perl scripts for each transactivation/repression domain model were run against the complete set of protein models from each genome ([Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online). 

Domain Co-occurrence Networks
-----------------------------

To visually represent the different domain architectures of the identified KLF/SP proteins, we created co-occurrence network maps.  For each species network map, all unique domain combination pairs were identified, summed, and divided by the total number of proteins (KLFs only, SPs only, or both KLF/SPs) to obtain a percentage of co-occurrence for a particular domain pair.  Composite network maps reflecting an organismal clade or grade were generated by the same procedure as for individual species maps.  The size of the circles (domains) reflects how often specific domains appear relative to the KLF-DBD, which has 100% representation.  Lines connecting two domains represent that unique combination pair.  Line weights represent the relative frequency of that specific pair combination.  A given domain pair combination was only counted once.  Network maps were visualized using Microsoft PowerPoint. 

Phylogenetic Analysis
---------------------

Putative filozoan KLF/SP and nonfilozoan amorphean KLF-DBD amino acid sequences were aligned using default settings of the MUSCLE alignment package in CLCBio ([@evv141-B33]).  To improve statistical support in our analyses, we concatenated the transactivation/repression domains with the KLF-DBD.  The highly variable 9aaTAD and LCR sequences were excluded from the alignment.  Individual aligned sequences ranged from 81 to 112 amino acids in length.  Each domain, including the KLF-DBD, can be distinguished as separate blocks within the alignment ([Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online).  For the nonfilozoan amorphean sequences, only the core three zinc fingers corresponding to the KLF-DBD model were included in the alignment.  Duplicate sequences were removed and are listed in [supplementary table S5](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online.  ProtTest v2.4 was used to determine the LG+I+G model as the best-fit model for protein evolution ([@evv141-B1]). 

Maximum likelihood (ML) analyses were performed using the MPI version of RAxMLv7.2.8 ([@evv141-B120]).  We executed 300 independent ML searches on randomized maximum parsimony starting trees using the standardized RAxML search algorithm, followed by comparison of likelihood values among all 300 resulting ML trees.  The final log-likelihood score of the best ML tree was −12974.002279.  One hundred bootstrap replicates were computed and applied to the best scoring ML tree.  ML bootstrap values are indicated on the ML tree ([supplementary fig.  S2](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online).  Bayesian analyses were performed with MrBayes3.2.5 ([@evv141-B95]).  We ran two independent 5 million generation runs of five chains each with default heating and with the "LG+I+G" amino acid model option.  The "average standard deviation of split frequencies" between the two runs was 0.048635.  This diagnostic is an indicator of how well the two runs converge.  A value below 0.01 is a strong indication of convergence, while a value between 0.01 and 0.05 is typically acceptable for convergence.  Additional convergence diagnostics were examined with AWTY ([@evv141-B79]), which was used to determine if a sufficient number of generations had been completed for posterior probabilities to stabilize, and to determine the amount of burn-in.  From 60,001 trees, 45,001 were sampled (25% burn-in was confirmed as adequate with Tracer v1.6; [@evv141-B91]) and used to create a consensus tree.  The runs reached stationarity, and adjusting the burn-in did not affect the topology of the tree.  Bayesian posterior probabilities (BPP) were calculated and are shown on the Bayesian consensus tree ([supplementary fig.  S3](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online).  FigTree v1.3.1 (<http://tree.bio.ed.ac.uk/software/figtree/>) was used for tree visualization. 

Results
=======

C2H2 Zinc Finger and KLF Identification
---------------------------------------

To better understand the evolution of the KLF/SP gene family and to gain insight into the genetic repertoire of transactivation/repression domains known to regulate disparate aspects of metabolism, and growth and development, we comprehensively searched for, identified, and characterized KLF/SP gene family complements from 48 eukaryotic genomes using a combination of hmms and custom perl models ([table 1](#evv141-T1){ref-type="table"}) with the assumption that the amorphean bikont split occurred at or near the origin of the Eukarya ([@evv141-B30]).  We found C2H2 zinc finger proteins highly represented in all 48 eukaryotic genomes.  The 81 amino acid KLF-DBD contains 3 highly conserved C2H2 zinc fingers separated by 2 highly conserved linker sequences ([supplementary fig.  S1*A*](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online).  We found this domain architecture to be restricted to the opisthokont lineage ([fig.  1](#evv141-F1){ref-type="fig"}).  Among a small number of nonfilozoans, represented by *Sphaeroforma* and several holomycotans, sequences containing the KLF-DBD motif were also found to possess additional zinc fingers.  Further analyses revealed the presence of a conserved aspartic acid residue in the second zinc finger of the KLF-DBD at position 44 (D~44~) among 397 putative KLF/SP filozoan sequences, the single exception being a ctenophore gene.  The D~44~ residue is critical for stabilizing proper KLF/SP DNA binding ([@evv141-B36]; [@evv141-B101]) and is notably absent from all nonfilozoan KLF-DBD motifs and the ctenophore sequence *MleKLFX* ([supplementary fig.  S1](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online).  Our phylogenetic analyses consistently recover a moderately supported clade that represents the 397 filozoan genes possessing the canonical D~44~ aspartic acid residue in the second zinc finger of the KLF-DBD to the exclusion of the remaining 13 nonfilozoan sequences, along with *MleKLFX* ([fig.  2](#evv141-F2){ref-type="fig"}).  We then searched the human genome with a refined KLF-DBD model including the canonical D~44~ residue.  We recovered the 26 human KLF/SPs plus Wilms Tumor 1 (WT1).  The WT1 gene contains an extra zinc finger and is phylogenetically distinct ([@evv141-B105]).  Therefore the high stringency D~44~ KLF-DBD model is specific and sensitive.  When the D~44~ KLF-DBD model is run against the eukaryotic genomes, we recover only putative KLF/SP protein sequences along with a small number of WT1 orthologs and WT1-like sequences that contain additional zinc fingers ([supplementary fig.  S1](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online).  F[ig]{.smallcaps}.  2.---Combined gene tree estimates for the concatenated KLF/SP data set using Bayesian criterion (MrBayes) and ML criterion (RAxML).  Gray node labels indicate congruent topology with BPP support = 84%.  Black node labels indicate congruent topology with BPP support ≥90%.  Clades collapsed to triangles indicate congruent topologies with BPP support ≥90%.  The single highly divergent ctenophore *MleKLFX* sequence clusters with nonfilozoan KLF-DBD presumably due to long-branch attraction.  Bayesian and ML trees with support values and branch lengths are available in [supplementary figs.  S2](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) and [S3](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online. 

Therefore, we operationally define bona fide KLF/SP orthologs as containing only 3 zinc fingers that conform to an 81 amino acid, 3 zinc finger KLF-DBD motif containing a canonical aspartic acid residue at position 44, D~44~ ([supplementary fig.  S1](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online).  Thus in our analyses KLF genes first appear in the filasterean *Capsaspora* while the SP subfamily first appears near the base of the metazoan lineage in sponges and is absent from ctenophores ([fig.  1](#evv141-F1){ref-type="fig"}).  The full set of KLF/SPs identified in this study with their corresponding accession numbers or protein/transcript identifiers are available in [supplementary table S1](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online.  Gene names given to each KLF/SP gene from this study can be found in [supplementary table S4](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online. 

### Bikonts

C2H2 zinc finger proteins were identified in both unicellular and multicellular bikont species.  However, the KLF-DBD was not present in any of the bikont genomes searched ([fig.  1](#evv141-F1){ref-type="fig"}). 

### Unicellular Amorpheans (Unikonts)

A phylogenetically diverse range of unicellular species, including amoebas and fungi, were used in this study.  Although C2H2 zinc finger proteins were identified in all unicellular amorpheans, the KLF-DBD was only found in *Fonticula* (sister to fungi), fungi (with the exceptions of *Mortierella* and *Encephalitozoon*), the ichthyosporean *Sphaeroforma*, the filasterean *Capsaspora*, and choanoflagellates ([fig.  1](#evv141-F1){ref-type="fig"}).  Among unicellular taxa bona fide KLF genes were only found in *Capsaspora* and choanoflagellates, two sister groups to the Metazoa ([fig.  1](#evv141-F1){ref-type="fig"}). 

### Metazoans

Within the metazoans C2H2 zinc finger proteins along with the KLF-DBD were present in all species.  However, neither representative ctenophore, *Mnemiopsis* or *Pleurobrachia*, show evidence for SP genes and both have fewer KLF genes as compared with other nonbilaterian metazoans.  Sponge, placozoan, and cnidarian KLF/SP complements are similar to bilaterian protostome lineages and early branching deuterostome lineages.  The total complement of KLF/SP genes substantially increases within the jawed vertebrates ([fig.  1](#evv141-F1){ref-type="fig"}) and teleosts have an average of 33 KLF/SPs, which is slightly higher than cartilaginous fishes (22) and tetrapods (25). 

Transactivation/Repression Domains
----------------------------------

The defining characteristic of the KLF/SP family, the highly conserved KLF-DBD, is located near the C-terminal region of nearly all identified KLF/SP proteins ([supplementary figs.  S1](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) and [S4](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online).  In contrast, N-terminal regions of KLF/SP proteins are generally less conserved.  These more variable N-terminal regions encode a variety of transactivation/repression domains, some of which have been well characterized in mammalian model systems ([@evv141-B123]; [@evv141-B70]).  We generated custom perl script models for the following transactivation/repression domains to facilitate their identification across filozoan KLF/SP proteins ([Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online): SID repressor domain, PVDLS repressor domain, R2/R3 repressor domains, the Btd box, and the SP box.  We screened for the presence of a conserved 9aaTAD ([@evv141-B85]).  We also identified several nonrandom LCRs implicated in transactivation/repression that are highly biased for a particular amino acid residue ([@evv141-B134]).  After determining the distribution of transactivation/repression domains within KLF/SPs ([fig.  3](#evv141-F3){ref-type="fig"}), we extended our search to include the full set of 48 eukaryotic genomes used in the study to better understand the representation of these domains in other proteins and to determine the extent domain shuffling may have played in the expansion of the KLF/SP gene family ([Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online).  F[ig]{.smallcaps}.  3.---Phylogenetic distribution of transactivation/repression domains and LCRs associated with KLF/SP proteins.  The**+**indicates the presence of the corresponding domain or LCR in at least one KLF/SP protein in the indicated taxa.  Only filozoan lineages containing bona fide KLF/SP proteins are shown.  An asterix indicates that RNA-seq data were used for that species.  Phylogeny is based on [@evv141-B31], [@evv141-B99], and [@evv141-B102]. 

### Low Complexity Regions

LCRs are nonrandom regions of protein sequences that are highly biased for a particular amino acid residue ([@evv141-B134]).  LCRs are commonly found in transcription factors ([@evv141-B35]), have been shown to influence transcriptional regulation ([@evv141-B38]), and are typically c-LCRs located centrally within the protein ([@evv141-B21]).  The composition of LCRs typically found in KLF/SP proteins includes serine/threonine (S/T)-rich, glutamine (Q)-rich, and proline (P)-rich regions.  S/T-rich regions are generally associated with enhanced transcriptional activation, whereas P-rich regions have been associated with transcriptional repression ([@evv141-B43]).  Q-rich regions, which are more frequently found among members of the SP subgroup, are known to interact with TAF~II~110 to activate transcription ([@evv141-B47]; [@evv141-B39]). 

We used four different algorithms (see Materials and Methods) to identify putative LCRs present in filozoan KLF/SP protein sequences ([supplementary table S2](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online).  We further required a putative LCR to be detected by a minimum of two methods for annotation as an LCR.  S/T-rich LCR regions occur most frequently and are found in at least one KLF/SP family member in all filozoan taxa except the poriferan *Oscarella* ([fig.  3](#evv141-F3){ref-type="fig"}).  P-rich LCR regions have lower representation among KLF/SPs, and are found in all filozoans except for *Ciona*, *Daphnia*, *Tribolium*, *Lottia*, *Trichoplax*, *Oscarella*, and *Monosiga* ([fig.  3](#evv141-F3){ref-type="fig"}).  KLF/SPs with Q-rich LCR regions were present in all filozoans except *Ciona*, *Tribolium*, *Capitella*, *Lottia*, cnidarians, *Oscarella*, ctenophores, and *Monosiga* ([fig.  3](#evv141-F3){ref-type="fig"}). 

### Nine-Amino-Acid Transactivation Domain

The 9aaTAD, first identified in yeast transcription factors ([@evv141-B85], is a short motif highly conserved throughout eukaryotes.  The 9aaTAD has been shown to interact with TAF9 of the RNA polymerase II holoenzyme and this domain motif has been identified in many transcription factors ([@evv141-B84]).  Using a 9aaTAD prediction tool ([@evv141-B85]), we found corresponding motifs in all filozoan KLF/SPs except *Capitella*, *Lottia*, and *Salpingoeca* ([fig.  3](#evv141-F3){ref-type="fig"}). 

### Buttonhead Box and SP Box

The Btd box and SP box are conserved domains that, in combination with the KLF-DBD, characterize the SP subfamily.  The Btd box was first identified in the *Drosophila* gene *buttonhead* ([@evv141-B132]) and the associated SP box was subsequently discovered in SP5 ([@evv141-B44]).  The Btd box motif is typically found just N-terminal of the KLF-DBD, while the SP box motif is located proximal to the N-terminus of the protein.  Although not definitive, there is evidence suggesting that the Btd box is involved in transactivation ([@evv141-B8]).  We examined filozoan KLF/SP sequences for the presence of the Btd and SP boxes.  Both Btd and SP box motifs are absent in unicellular KLFs, absent in both ctenophore KLFs, but present in all other metazoan phyla ([fig.  3](#evv141-F3){ref-type="fig"}).  Uniquely within the poriferans, the SP box was not found in any *Amphimedon* KLF/SP, but was identified in *Oscarella* KLF/SPs. 

Using the same domain models to screen the complete set of 48 eukaryotic whole genomes, we identified Btd box motifs in a number of genes other than the KLF/SP family with the notable exception of the fungi, *Encephalitozoon* and *Saccharomyces* ([Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online).  In contrast, the distribution of the SP box motif is highly restricted.  This domain motif first appears coincident with the SP subfamily in poriferans.  The SP box is present in one out of two SP genes in *Oscarella* but is not detected in the two SP genes identified in the *Amphimedon* genome.  Excluding the ctenophores, metazoan genomes outside the poriferans have SP box domain motifs among a small set of genes to the exclusion of the SPs ([fig.  6](#evv141-F6){ref-type="fig"}). 

### Sin3a Interacting Domain

The Sin3a protein acts as a transcriptional repressor, and is able to recruit and bind histone deacetylases ([@evv141-B62]; [@evv141-B107]).  Mammalian KLF9, 10, 11, 13, 14, and 16 are known to interact with the Sin3a protein through a Sin3a interacting domain (SID) which binds the PAH domain of the Sin3a protein ([@evv141-B49]; [@evv141-B13]; [@evv141-B23]; [@evv141-B115]; [@evv141-B55]).  No SID motif was detected in unicellular KLF/SP sequences.  Within metazoans a SID-containing KLF was identified in all taxa except ctenophores and the protostomes *Capitella*, *Lottia*, and *Drosophila* ([fig.  3](#evv141-F3){ref-type="fig"}).  To gain further insight into the evolutionary history of the SID motif, we searched the complete set of 48 whole genomes and identified a conserved SID in a number of genes in all representative eukaryote genomes ([Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online). 

### R2 and R3 Domains

KLF10 and KLF11 proteins contain R2 and R3 repressor domains in combination with an SID motif ([@evv141-B22]).  We searched filozoan KLF/SP complements and were only able to identify R2 and R3 domains in representative vertebrate taxa ([fig.  3](#evv141-F3){ref-type="fig"}).  We then searched the complete set of 48 whole genomes for the presence of the R2 and R3 domains.  Our search revealed the presence of an R3 domain in all representative eukaryotic genomes, while the R2 domain was restricted to vertebrate KLF10/11 genes ([Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online). 

### PVDLS Domain

KLF3, KLF8, and KLF12 can corepress transcription through interaction with the C-terminal binding protein mediated by the PVDLS domain ([@evv141-B24]; [@evv141-B127]; [@evv141-B50]; [@evv141-B129]).  A PVDLS domain was identified in at least one KLF/SP in all jawed vertebrates and in *Drosophila* ([fig.  3](#evv141-F3){ref-type="fig"}).  To help resolve this relationship between the fly PVDLS containing KLF and the vertebrate KLF3/8/12 genes, we searched additional protostome genomes for the presence of KLF genes containing the PVDLS motif.  Curiously, the association of this motif with KLF genes was not detected in other Drosophilid species, but was identified in two hymenopterans *Apis millifera* and *Nasonia vitripennis* (data not shown).  Analysis of the complete set of 48 whole genomes revealed the PVDLS domain in multiple genes in all representative eukaryote species ([Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online). 

Co-occurrence Networks
----------------------

Co-occurrence networks are visual representations of domain pair occurrences within a given protein or protein family.  Typically, domains are represented as circles, and lines connect domains that appear together within a protein.  These maps can be useful for visualizing the frequencies of certain domains occurring with each other, that is, domain architectures.  To explore differences in domain architecture complexity, we generated domain co-occurrence maps for different KLF/SP domain networks ([fig.  4](#evv141-F4){ref-type="fig"}).  These maps can also show the general N-terminal to C-terminal relationship between different domains ([fig.  4](#evv141-F4){ref-type="fig"}*A*).  Our co-occurrence maps indicate the frequency of unique domain pairs as well as how often an individual domain appears within a given network in extant species ([fig.  4](#evv141-F4){ref-type="fig"}*B*--*I*).  Composite networks of larger taxonomic groupings represent the consensus map of all included species.  The composite networks are additive and thus do not compensate for missing data due to poor genome annotations.  Importantly, these co-occurrence maps are not intended to be ancestral reconstructions.  They are used to show observed domain relationships in, or among groupings of, extant species.  Nonetheless, evolutionary inferences can be drawn from these network maps when combined with inferences of ancestral presence or absence of discrete transactivation/repression domains as represented in [figure 6](#evv141-F6){ref-type="fig"}.  F[ig]{.smallcaps}.  4.---KLF/SP protein domain co-occurrence networks.  In all networks, each circle represents a transactivation/repression domain or an LCR.  A line connecting two domains indicates a co-occurrence of those two domains.  Domains are arranged in approximately the same 5′--3′ spatial orientation as they appear encoded in KLF/SP sequences. (* A*) General network diagram showing connectivity and unidirectional spatial relationships between transactivation domains among filozoan KLF/SPs.  Blue arrows represent connectivity upstream of the KLF-DBD; the gold arrow represents connectivity downstream of the KLF-DBD. (* B--I*) KLF/SP co-occurrence networks from different taxonomic groups.  Circle size indicates the relative frequency of occurrence in the network, with the KLF-DBD always representing 100%.  Circle color follows the same convention as seen in [figure 3](#evv141-F3){ref-type="fig"}.  Repeated domains were counted as occurring only once.  Lines connecting circles indicate the presence of that specific domain pair co-occurrence in at least one KLF/SP.  Line width indicates the frequency of domain pair co-occurrence.  Only LCR domains which are found N-terminal of the KLF-DBD are represented in these networks ([supplementary fig.  S4](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online). (* B*) Complete filozoan KLF/SP network. (* C*) Representative unicellular KLF/SP network. (* D*) KLF/SP network from nonbilaterian metazoans. (* E*) Invertebrate bilaterian KLF/SP network. (* F*) Vertebrate KLF/SP network. (* G-H*) Representative ctenophoran and poriferan KLF/SP networks for comparison with each other and with the network in *D*. (* I*) *Ciona* KLF/SP network for comparison with the networks in *E* and *F*. (* J*, *K*) Co-occurrence network maps for the KLF subfamily and SP subfamily mapped onto the filozoan phylogeny ([@evv141-B31]; [@evv141-B99]) for evolutionary comparison.  Each network represents a composite for the taxonomic group indicated. (* J*) Co-occurrence maps for domains found in the KLF subfamily. (* K*) Co-occurrence maps for domains found in the SP subfamily.  The unicellular filozoan genomes and ctenophore genomes do not contain SP genes. 

The unicellular KLF/SP network (*Capsaspora* + choanoflagellates) is one of the least complex with only LCR domains linked to the KLF-DBD ([fig.  4](#evv141-F4){ref-type="fig"}*C* ).  Within metazoans, there appears to be a gradient of network complexity; the nonbilaterian network ([fig.  4](#evv141-F4){ref-type="fig"}*D*) is less complex than the invertebrate bilaterian network ([fig.  4](#evv141-F4){ref-type="fig"}*E*), and both of these are less complex than the vertebrate network ([fig.  4](#evv141-F4){ref-type="fig"}*F*).  A small number of lineage-specific networks showed a significant departure from larger composite networks.  For example, the Ctenophoran network ([fig.  4](#evv141-F4){ref-type="fig"}*G*), in contrast to both the Poriferan network ([fig.  4](#evv141-F4){ref-type="fig"}*H*) and larger nonbilaterian composite network ([fig.  4](#evv141-F4){ref-type="fig"}*D*), is composed of only S/T-rich and P-rich LCRs, with the P-rich LCRs + KLF-DBD domain pair occurring at greater frequency.  The Poriferan network ([fig.  4](#evv141-F4){ref-type="fig"}*H*) bears substantially greater similarity to the more inclusive nonbilaterian composite network ([fig.  4](#evv141-F4){ref-type="fig"}*D*).  The urochordate *Ciona* also presents an interesting departure from composite networks ([fig.  4](#evv141-F4){ref-type="fig"}*I*), diverging dramatically from the invertebrate bilaterians in having S/T-rich LCRs as the most prevalent domain linked to the KLF-DBD.  The relationship of Ctenophora and Porifera to other metazoans shown in [figure 4](#evv141-F4){ref-type="fig"}*J* and *K* are based on recent hypotheses of metazoan phylogeny ([@evv141-B31]; [@evv141-B99]). 

KLF Domain Architecture and Phylogenetics
-----------------------------------------

All KLF/SP gene family members share homology at the conserved KLF-DBD, typically located toward the C-terminus.  The N-terminal regions share very little similarity across the entire family.  Distinct KLF/SP subgroups can, however, be defined based on comparable structure and function ([@evv141-B70]).  These subgroups share domain architectures defined by unique combinations of transactivation/repression domains occurring with the KLF-DBD ([fig.  5](#evv141-F5){ref-type="fig"}; [supplementary table S3](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online).  Domain architectures have been described for mammalian KLF/SP genes ([@evv141-B54]; [@evv141-B123]; [@evv141-B70]; [@evv141-B7]).  Historically, KLFs have been divided into five domain architectures and named according to the genes in each class: KLF1/2/4, KLF6/7, KLF3/8/12, KLF9/13, and KLF10/11.  Similarly, the SPs have been divided into two architecture classes: SP1--4 and SP5--9.  Using this existing classification scheme as a foundation, we determined the distribution of domain architectures for all filozoan KLF/SP genes.  The 9aaTAD was not used in determining domain architecture as it does not contribute to any unique domain combination.  F[ig]{.smallcaps}.  5.---Phylogenetic distribution of explicit domain architectures represented among KLF/SP proteins.  The key at lower left identifies LCRs and transactivation/repression domains used to determine domain architectures.  The protein schematics along lower right represent the particular combinations of domains and LCRs with the KLF-DBD that define each specific KLF/SP protein architecture.  All groups, except for the ancient unicellular KLF architecture recovered, are named according to established human KLF/SP paralogy groups that conform to each specific architecture.  The three C-terminal zinc fingers of the KLF-DBD are indicated with grey boxes labeled zf1, zf2, and zf3.  Architecture schematics are not to scale.  F[ig]{.smallcaps}.  6.---Inferred relationships between key events during the evolution and expansion of the KLF/SP gene family.  Symbol key is at upper left.  Colored rectangles represent the origin of particular transactivation/repression domains or LCRs co-occurring with the KLF-DBD ([fig.  4](#evv141-F4){ref-type="fig"}).  Yellow hexagons represent the origin of specific KLF/SP domain architectures ([fig.  5](#evv141-F5){ref-type="fig"}).  A black X over a hexagon represents the loss of specific domain architecture.  Colored triangles represent the presence of specific transactivation domain motifs within whole eukaryote genomes to the exclusion of the KLF/SP gene family. (* A*) We infer the origin of the KLF-DBD in the opisthokont stem lineage prior to the divergence of the Holomycota.  However, bona fide KLF gene architectures do not appear until the divergence of the filozoan lineage (KLF origin).  The ancient unicellular KLF domain architecture is not recovered in metazoan lineages.  The ancient PVDLS, SID, Btd box, and R3 domains were recovered, to the exclusion of KLF/SPs, in all eukaryote genomes searched.  Notably, the Btd box was not recovered in *Saccharomyces* and *Encephalitozoon* fungal genomes.  Our analysis suggests that the origin of the SP subfamily is in the metazoan stem lineage prior to the divergence of the poriferans; it is not present in the ctenophorans.  The SP box motif only appears in SP genes in poriferans and is not found in additional genes until the divergence of *Trichoplax*.  The R2 repressor domain appears to be a de novo innovation restricted to KLF genes in the vertebrate stem lineage, contributing to the KLF10/11 architecture class.  Composite domain co-occurrence maps for each taxonomic group are shown to the right of the tree.  Representative examples of putative domain shuffling events during the evolution and expansion of the KLF/SP gene family. (* B*) An ancient Btd box and a metazoan SP gene may have contributed to the origin of the SP gene subfamily early in metazoan evolution. (* C*) An ancient SID likely combined with a pre-existing ancestral KLF gene to form the KLF9/13 group, also early in metazoan evolution. (* D*) An ancient PVDLS domain combined with a pre-existing ancestral KLF gene to form the KLF3/8/12 group.  We infer an independent convergent acquisition of the PVDLS domain within a KLF gene in the Protostomia lineage (see Discussion).  Domain icon colors are the same as [figure 5](#evv141-F5){ref-type="fig"}. 

We performed ML and Bayesian phylogenetic analyses on an alignment of nonfilozoan amorphean and filozoan concatenated sequences that included the KLF-DBD and transactivation domains ([supplementary figs.  S2](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) and [S3](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online).  The overall topologies of the ML and Bayesian trees share significant overlap, albeit with only moderate support for deeper nodes, for example, BPP = 0.84 ([fig.  2](#evv141-F2){ref-type="fig"}).  A number of informative clades with congruence across both ML and Bayesian analyses were recovered.  We primarily focused on overlapping nodes with BPP ≥ 0.90.  Most of these nodes are consistent with explicit domain architectures and well-accepted inferences of metazoan phylogeny.  Many of the well-supported orthologous sequence clades are composed of vertebrates with corresponding nonvertebrate putative orthologs diverging prior to the highly supported vertebrate nodes ([fig.  2](#evv141-F2){ref-type="fig"}; [supplementary figs.  S2](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) and [S3](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online). 

Among the unicellular KLF representatives, the *Capsaspora* and *Salpingoeca* domain architectures consist of alternating S/T- and Q-rich LCRs ([fig.  5](#evv141-F5){ref-type="fig"}) with scattered P-rich LCRs at low frequency ([fig.  4](#evv141-F4){ref-type="fig"}*C*).  Although the S/T- and Q-rich LCR organization superficially resembles the architecture found among SP1--4 genes, these unicellular KLFs lack the characteristic Btd box and SP box motifs that define the SPs.  Therefore, we created a unicellular specific KLF architecture class ([fig.  5](#evv141-F5){ref-type="fig"}).  Although a single unicellular KLF clade composed of *CowKLFa*, *SroKLFa*, and *MbrKLFa* branch sister to a highly supported cnidarian + bilaterian KLF15 clade (BPP = 0.99), the remaining unicellular KLF sequences are divergent and characterized by long branches ([fig.  2](#evv141-F2){ref-type="fig"}; [supplementary figs.  S2](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) and [S3](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online). 

The KLF6/7 and KLF1/2/4 classes are defined by the presence of mostly S/T-rich or P-rich LCRs, respectively.  KLF6/7 domain architectures were identified in *Monosiga* and all metazoans except for *Oscarella* and *Lottia* ([fig.  5](#evv141-F5){ref-type="fig"}).  The single *Monosiga* KLF gene superficially consists of a KLF6/7 domain architecture; however, phylogenetically this sequence groups with other unicellular KLFs. * Amphimedon* sequences with KLF6/7 architecture grouped more closely with KLF1/2/4 and KLF5, albeit with low statistical support.  Our phylogenetic analysis recovered a well-supported clade (BPP = 0.91) that includes bilaterian KLF6/7 along with representative cnidarian and *Trichoplax* KLF6/7.  The KLF1/2/4 class was identified in all metazoans except *Oscarella*, *Trichoplax*, *Lottia*, *Drosophila*, *Daphnia*, and *Ciona* ([fig.  5](#evv141-F5){ref-type="fig"}).  Our phylogenetic analyses recovered a KLF 1/2/4 clade including invertebrate bilaterians (BPP = 0.96), a teleost-specific KLF1 clade (BPP = 1), and a mammalian KLF1 clade (BPP = 1; [fig.  2](#evv141-F2){ref-type="fig"}). 

The KLF9/13 class possesses a SID near the N-terminal region and typically contains one or more c-LCRs located between the SID and KLF-DBD.  The related KLF10/11 class shares the SID with KLF9/13, but uniquely contains the R2 and R3 repressor domains situated between the SID and the KLF-DBD.  The R2 and R3 repressor domains are generally flanked by one or more P-rich or S/T rich c-LCRs.  Based on phylogenetic analyses, all metazoans, with the notable exception of ctenophores, have at least one KLF9/13 gene ([supplementary figs.  S2--S4](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online).  Mammalian KLF14 and KLF16 form a unique, highly supported (BPP = 1) clade within the 9/13 group supporting the assignment of these orthologs to a mammalian-specific lineage ([fig.  2](#evv141-F2){ref-type="fig"}; [supplementary fig.  S5](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online).  The KLF10/11 class has a more restricted distribution and is exclusively found in vertebrates ([fig.  5](#evv141-F5){ref-type="fig"}).  All SID containing KLF sequences form a clade with low support (BPP = 0.75) ([supplementary figs.  S2](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) and [S3](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online). 

Finally, the KLF3/8/12 class is characterized by a PVDLS domain and often has additional c-LCRs between the PVDLS and the KLF-DBD.  According to domain architecture, this class is found in all the jawed vertebrates and in the Endopterygota as represented in this study by *Drosophila* ([fig.  5](#evv141-F5){ref-type="fig"}).  Our phylogenetic analyses recover a KLF3/8/12 clade with high support (BPP = 0.95; [fig.  2](#evv141-F2){ref-type="fig"}).  Notably, the *Drosophila* KLF sequence that contains a PVDLS domain (*CG42741*) falls outside of the KLF3/8/12 clade which consists exclusively of vertebrate sequences.  Within this clade the KLF8 orthologs are recovered with BPP = 1.  Despite poor annotation in the lamprey genome, we identified a lamprey sequence within the vertebrate KLF3/8/12 clade ([supplementary figs.  S2](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) and [S3](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online). 

The SP subfamily is defined by the presence of a Btd box 5′ proximal to the KLF-DBD and generally possesses an additional SP box located near the N-terminus.  The SP factors can be further separated into two domain architecture classes, SP1--4 and SP5--9 ([@evv141-B14]; [@evv141-B7]).  The SP1--4 class typically has S/T-rich c-LCRs adjacent to Q-rich c-LCRs.  The SP5--9 class is characterized by having either predominately S/T, or with lesser frequency, P-rich c-LCRs.  Notably, the SP subfamily is absent in both ctenophore genomes.  The SP1--4 domain architecture is found in all metazoans except *Oscarella*, cnidarians, *Capitella*, *Lottia*, *Tribolium*, and *Ciona*.  The SP5--9 domain architecture is found in all metazoan genomes excluding ctenophores, *Drosophila*, and *Daphnia* ([fig.  5](#evv141-F5){ref-type="fig"}).  Phylogenetic analyses recover moderate support (BPP = 0.84) for a metazoan SP clade ([fig.  2](#evv141-F2){ref-type="fig"}; [supplementary figs.  S2](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) and [S3](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online).  Well-supported SP clades (BPP ≥ 0.9) include vertebrate SP2, SP6, and SP7. 

Discussion
==========

KLF/SP Gene Family Origins in Eukarya
-------------------------------------

The evolution of the KLF/SP family and associated transactivation/repression domains is represented in [figure 6](#evv141-F6){ref-type="fig"}.  Although C2H2 zinc finger domains are ubiquitous across a wide range of eukaryotes ([@evv141-B28]), the KLF-DBD first appears in the Holomycota ([figs.  1](#evv141-F1){ref-type="fig"} and [6](#evv141-F6){ref-type="fig"}).  Historically, the single criteria for defining the KLF/SP gene family has been the presence of a highly conserved KLF-DBD composed of three C2H2 zinc fingers each separated by a seven amino acid linker region ([@evv141-B70]).  The first 2 zinc fingers are 23 amino acids long (from C to H), while the third zinc finger is only 21 amino acids long.  Our analyses reveal that all filozoan KLF/SP sequences recovered except one, the highly divergent *MleKLFX* ctenophore sequence, contains a canonical aspartic acid residue in the second zinc finger, D~44~, which is critical for proper DNA binding ([supplementary fig.  S1](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online) ([@evv141-B36]; [@evv141-B101]).  All nonfilozoan KLF-DBD containing genes lacked this diagnostic residue and possess additional zinc fingers.  Furthermore, none of the transactivation/repression domains considered in this study were found in the nonfilozoan KLF-DBD containing sequences, despite several relevant domains being well represented in all eukaryotic genomes examined ([fig.  6](#evv141-F6){ref-type="fig"}, [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online).  Thus the nonfilozoan KLF-DBD sequences were not classified as bona fide KLF/SP orthologs.  Moreover, our phylogenetic analyses recovered a moderately supported filozoan clade that possesses the canonical D~44~ residue to the exclusion of all nonfilozoan sequences and the highly divergent *MleKLFX* ctenophore sequence ([fig.  2](#evv141-F2){ref-type="fig"}).  We operationally define a bona fide KLF/SP gene as containing only 3 zinc fingers conforming to 81 amino acid, 3 zinc finger KLF-DBD motif typically found in the C-terminal region of the parent sequence and containing a canonical aspartic acid residue at position 44, D~44~, following the general consensus sequence C-X~4~-C-X~12~-H-X~3~-H-X~7~-C-X~4~-C-X~7~-**D**-X~4~-H-X~3~-H-X~7~-C-X~2~-C-X~12~-H-X~3~-H.

Therefore, the KLF-DBD likely has its origins in the opisthokont stem lineage prior to the divergence of the Holomycota ([figs.  1](#evv141-F1){ref-type="fig"} and [6](#evv141-F6){ref-type="fig"}).  In our study, KLF genes first appear in the filasterian *Capsaspora*, while the SP subfamily is restricted to metazoans, notably excluding the ctenophores.  Thus we infer the origin of KLF genes in the filozoan stem lineage and a later origin of the SP genes after the divergence of the metazoan lineage ([figs.  1](#evv141-F1){ref-type="fig"} and [6](#evv141-F6){ref-type="fig"}).  Our phylogeny assumes an early divergence of the ctenophore lineage from the metazoan stem ([@evv141-B31]; [@evv141-B99]).  An alternate view of early animal phylogeny ([@evv141-B83]; [@evv141-B77]) would infer the origin of SP genes prior to the divergence of the poriferans and assume a subsequent loss of SP orthologs in the ctenophore lineage.  Under either scenario of early metazoan lineage divergence, the appearance of KLF genes precedes the origin of the SP gene subfamily. 

Apart from the KLF-DBD, the most common sequence features are the presence of one or more S/T-, P-, or Q-rich c-LCRs ([figs.  3--5](#evv141-F3 evv141-F4 evv141-F5){ref-type="fig"}; [supplementary fig.  S4](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online).  Our results highlight promiscuous variation in LCR composition and length between KLF/SP gene family members ([supplementary fig.  S4](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online).  In contrast to earlier studies with more restricted sampling, we are unable to discriminate phylogenetic relationships between KLF/SP architectural classes composed of explicit c-LCR compositional types linked to the KLF-DBD.  For example, many representative genes nested within the KLF1/2/4 and KLF6/7 clades include all three predominant c-LCR compositional types ([figs.  2](#evv141-F2){ref-type="fig"} and [5](#evv141-F5){ref-type="fig"}; [supplementary fig.  S4](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online).  Variability, expansion, and extinction of LCRs have been associated with gene conversion due to mismatch repair of DNA heteroduplexes ([@evv141-B90]) and higher rates of recombination relative to flanking sequences due to unequal crossing over and replication slippage ([@evv141-B139]).  Both of these core genetic mechanisms would contribute to our observations of significant LCR homoplasy within and between orthologous domain architectures. 

Our analyses of the frequency of appearance of explicit compositionally biased c-LCRs and their frequency of co-occurrence with flanking transactivation domains suggest functional consequences.  For example, S/T- and Q-rich LCRs are typically acidic and associated with transcriptional activation, whereas P-rich regions, with their bulky side chains, are typically neutral and associated with context-dependent repression.  The high frequency of S/T-rich c-LCRs suggests an ancestral KLF that was likely involved in transcriptional activation.  Despite high rates of gene conversion and unequal crossover associated with LCRs, we observe P-rich c-LCRs at low frequency.  Conversely, among the KLF/SP genes that retain P-rich c-LCRs, our analyses also show an increase in the frequency of connectivity with other repressor motifs in metazoans ([fig.  4](#evv141-F4){ref-type="fig"}).  This suggests the evolution of coordinated multidomain repression during the expansion of the KLF/SP gene family.  For example, known repressor domains such as the SID and PVDLS motifs are typically located 5′ of c-LCRs and show a similar pattern of increased frequency during metazoan lineage diversification ([figs.  4](#evv141-F4){ref-type="fig"} and [6](#evv141-F6){ref-type="fig"}).  This pattern suggests that within the KLF/SP gene family, despite a selective preference for acidic c-LCRs as evidenced by their high frequency, the evolution of a repertoire of transcriptional repressor combinations through the differential pairing of c-LCRs with repressor motifs occurred.  Our analyses further suggest that these shuffling events may be associated with cell type diversification in metazoans.  For example, chordates possess \>100 distinct cell types, whereas earlier diverging metazoan lineages range from 4 to 59 distinct cell types ([@evv141-B16]).  This pattern mirrors the increased frequency of repressor domain connectivity observed across members of the KLF/SP gene family in metazoans ([figs.  4--6](#evv141-F4 evv141-F5 evv141-F6){ref-type="fig"}). 

The total complement of KLF/SP genes across nonbilaterian metazoans, protostomes, and invertebrate deuterostomes shows only a ∼2-fold variance of 6--12 genes.  Vertebrates have an average number of 24 KLF/SPs, likely due to the 2 rounds of whole genome duplication (WGD) in the vertebrate stem lineage ([supplementary fig.  S5](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online).  Despite the additional round of WGD in the teleosts, their KLF/SP gene complement is not substantially greater than other jawed vertebrates.  It seems that KLF/SP gene duplicates were not necessarily preferentially retained, which has been shown to be the case with other transcription factors ([@evv141-B28]).  The lamprey *Petromyzon marinus* has an atypically depauperate KLF/SP complement as compared with other vertebrates.  It has been shown that during embryogenesis lampreys can lose ∼20% of germline DNA in somatic tissues due to genomic rearrangements ([@evv141-B111], [@evv141-B112]).  The lamprey genome was derived completely from somatic tissue, thus the true *P.marinus* KLF/SP complement may be underrepresented in our analyses ([@evv141-B113]). 

Ctenophores have a reduced KLF/SP complement compared with other nonbilaterian metazoans and lack members of the SP subfamily altogether ([fig.  1](#evv141-F1){ref-type="fig"}).  This pattern of gene underrepresentation in ctenophores relative to other metazoans has been observed in a number of studies ([@evv141-B98]; [@evv141-B68]; [@evv141-B74]; [@evv141-B102]).  Among the filozoan KLF orthologs, a single member, the highly divergent ctenophore gene *MleKLFX*, possesses atypical domain architecture in which a KLF-DBD lacking the canonical D~44~ residue is located in the N-terminal region of the protein instead of the C-terminal region ([supplementary fig.  S4](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online).  This gene consistently groups with the nonfilozoan amorphean KLF-DBD sequences presumably due to the effects of long-branch attraction ([fig.  2](#evv141-F2){ref-type="fig"}; [supplementary figs.  S2](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) and [S3](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online). 

The Btd box is found in SPs and along with the KLF-DBD represents the defining characteristic of the SP subfamily.  The Btd box was initially described in the *Drosophila btd* gene which has high similarity to mouse SP1 and SP3 ([@evv141-B132]).  Most of the amino acid sequence conservation resides in the zinc fingers, and *btd* is often classified as an SP homolog in flies.  However, the second zinc finger of the *btd* gene is slightly different than the corresponding zinc finger in the KLF-DBD.  The second zinc finger of *btd* is only 21 amino acids long, lacking 2 residues between the 2 C residues.  According to this structural difference, *btd* would be classified in the Zif268/EGR-1 family of zinc finger transcription factors ([@evv141-B51]).  However, phylogenetic analysis of *btd* shows that it is more similar to SPs than to EGR-1 family genes and *btd*-like orthologs were not found in other organisms (data not shown).  These suggest that at some point along the stem lineage leading to *Drosophilia*, there was a lineage-specific deletion of 2 amino acid residues in the second zinc finger of *btd*.  Notably, the loss of these two residues does not hinder *btd* from binding transcription factor binding sites similar to SP factors ([@evv141-B132]). 

Recently, three SPs in the cnidarian *Nematostella* and the placozoan *Trichoplax* were identified and grouped into three separate clades: SP1--4, SP5/*btd*, and SP6-9 ([@evv141-B100]) leading to a proposed ancestral complement of three SPs in the Metazoa.  In our study, we recovered two SPs in sponges, two SP genes in *Trichoplax*, and four SP genes in *Nematostella*.  The *Trichoplax btd* gene and SP6--9 gene from [@evv141-B100] match with two *Trichoplax* SPs from our study ([supplementary table S1 and](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) [supplementary fig.  S4](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online).  However, the *Trichoplax* SP1--4 gene from [@evv141-B100] matches with the *Trichoplax* KLF9/13 gene in our study.  We were able to extend this *Trichoplax* sequence and identified a SID motif in the previously unannotated 5′ end of this gene model thus confirming the placement of this sequence within the KLF9/13 group ([supplementary fig.  S4](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online).  Therefore, our analyses suggest an ancestral complement of a single SP1--4 ortholog and a single SP5--9 ortholog in the early metazoan stem lineage followed by an expansion after the divergence of the Placozoa ([fig.  1](#evv141-F1){ref-type="fig"}). 

Our analyses also suggest that two KLF9/13 genes were present early in Metazoa ([figs.  5](#evv141-F5){ref-type="fig"} and [6](#evv141-F6){ref-type="fig"}; [supplementary figs.  S2](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) and [S3](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online).  In mammals, two additional genes with similar domain architecture, KLF14 and KLF16, are present and a recent study provided evidence that KLF14 evolved from KLF16 through a retrotransposon event (RTE) within the mammalian stem lineage ([@evv141-B81]).  Thus, we infer that KLF16 arose by a tandem duplication of one of the two ancestral KLF9/13 genes in the mammalian stem lineage followed by KLF14 evolution by RTE of KLF16 ([supplementary fig.  S5](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online). 

Transactivation/Repression Domains Show Unique Evolutionary History
-------------------------------------------------------------------

Our results highlight that individual transactivation/repression domains associated with KLF/SP transcription factors have unique evolutionary histories ([fig.  6](#evv141-F6){ref-type="fig"}).  The Btd box, PVDLS, R3, and SID are ancient domain motifs present in other genes in all representative eukaryote genomes in this study.  The notable exceptions being the *Saccharomyces* and *Encephalitozoon* fungal genomes in which no Btd box domains were detected.  These four ancient domains first appear in differential combinations with the KLF-DBD in KLF/SP genes after the divergence of the metazoans and correlate strongly with the observed expansion of domain architecture repertoires ([fig.  5](#evv141-F5){ref-type="fig"}).  Given the phylogenetic distribution of domain architectures associated with the putative acquisition of ancient motifs by domain shuffling, we infer that particular domain architectures stem from unique independent ancestral shuffling events along the metazoan stem lineage ([figs.  5](#evv141-F5){ref-type="fig"} and [6](#evv141-F6){ref-type="fig"}). 

Our analyses also infer two instances of de novo domain origin within the KLF/SP gene family ([fig.  6](#evv141-F6){ref-type="fig"}).  The SP subfamily appears early in metazoan evolution in the Poriferan lineage ([fig.  1](#evv141-F1){ref-type="fig"}) coincident with the appearance of the SP box motif ([fig.  3](#evv141-F3){ref-type="fig"}).  However, in contrast to the Btd box motif, the SP box motif has a very limited genomic and phylogenetic distribution.  In sponges, the SP box motif is uniquely associated with the SP genes.  In later diverging lineages of metazoans, the SP box is associated with a small number of genes in addition to the KLF/SP family.  Thus the SP genes provide an example of a metazoan-specific multidomain protein that consists of both ancient domains, including the Btd box and the KLF-DBD, coupled with the de novo origin of a metazoan-specific domain motif, the SP box.  Another example of de novo domain motif origin within the KLF/SP family is the vertebrate-specific KLF10/11 genes ([fig.  6](#evv141-F6){ref-type="fig"}; [supplementary fig.  S5](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online).  This architectural class is composed of ancient SID, R3, and KLF-DBD domains combined with a vertebrate-specific R2 domain.  Our exhaustive search uncovered no R2 domains in any invertebrate genome or outside of the KLF10/11 genes within the vertebrate genomes in this study.  Interestingly, the R3 domain is present in all eukaryote genomes but there is no evidence for it being shuffled into KLFs until much later than the SID motif ([fig.  6](#evv141-F6){ref-type="fig"}). 

The PVDLS domain represents an intriguing case of putative convergence.  In our analyses, the PVDLS domain appears in vertebrate KLFs defining the 3/8/12 architecture group and in a lone *Drosophila* KLF (*CG42741*).  Our phylogenetic analyses suggest, however, that the fly sequence is not nested within the highly supported KLF3/8/12 clade ([supplementary figs.  S2](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) and [S3](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1), [Supplementary Material](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) online).  To help elucidate the incongruous relationship between the representative fly gene containing a PVDLS motif and the PVDLS containing vertebrate KLF3/8/12 class, we searched a number of additional protostome genomes for the presence of KLF genes containing the PVDLS motif (data not shown).  Our search yielded only two other instances, both within the hexapods, of a KLF gene also containing a PVDLS motif: The hymenopterans *A. mellifera* and *N.vitripennis*.  Exhaustive searches of 12 other Drosophilid species did not uncover any PVDLS containing KLFs.  Based on these additional results, we infer that a putative ancestral KLF3/8/12 gene most likely evolved early in the vertebrate stem lineage, whereas the PVDLS motif was likely also convergently acquired in the hexapod lineage leading to the Endopterygota. 

Conclusions
===========

Our analysis across 48 eukaryotic genomes illuminates the origin and evolutionary history of the KLF/SP gene family.  We also identify and characterize associated transactivation/repression domains, including LCRs, enabling us to develop models of KLF/SP domain co-occurrence evolution.  By extending our domain search to include entire proteomes, we find evidence for a complex intersection of domain shuffling, gene duplication, and de novo domain evolution as the primary mechanisms for the diversification of the KLF/SP gene family across the Metazoa.  Our results uncover a pattern of an increased frequency of repressive domain connectivity repertoires (P-rich LCRs, SID, R2/R3, and PVDLS domains) in the KLF/SP gene family among metazoans suggesting a role in mediating diverse transcriptional repression activity.  Our phylogenetic results further suggest that the expansion of the KLF/SP gene family mirrors increased cell type diversity during metazoan lineage diversification.  The expansion and diversification of the KLF/SP gene family within the Metazoa may thus reflect the accumulation of differential transcriptional repression strategies associated with the development of extensive repertoires of cell types required to support complex tissues. 

Supplementary Material
======================

[Supplementary figs.  S1--S5](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) and [tables S1--S5](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) are available at Genome Biology and Evolution online (<http://www.gbe.oxfordjournals.org/>).  Additional online [supplementary materials](http://gbe.oxfordjournals.org/lookup/suppl/doi:10.1093/gbe/evv141/-/DC1) can be viewed at <https://goo.gl/dbdBil>. 

###### Supplementary Data

This work was supported by startup funds from the University of Miami College of Arts and Sciences to W.E.B. J.S.P. was supported by the University of Miami College of Arts and Sciences.  C.E.S. was supported by the Intramural Research Program of the National Human Genome Research Institute, National Institutes of Health.  The authors thank James Baker, Allie Graham, Isaac Skromne, Lauren Vandepas, members of the Browne laboratory, and two anonymous reviewers for their thoughtful comments and critical reading of the manuscript. 

[^1]: **Associate editor:** Mar Alba

[^2]: **Data deposition:** Sequences associated with this project have been deposited at GenBank under accession numbers KJ576841 and KJ576842.