Genome-wide identification of Tribolium dorsoventral patterning genes

In the fruit fly, Drosophila melanogaster, Toll signalling is crucial for both dorsoventral (DV) axis formation and innate immunity (Leulier and Lemaitre, 2008; Stein and Stevens, 2014). Although the immune function of Toll is conserved from flies to vertebrates, a DV patterning function of Toll has not been observed outside the insects. Indeed, analysis of DV patterning in more basally branching insects indicates that during insect evolution Toll signalling progressively gained the complex morphogen function found in modern flies (Lynch and Roth, 2011; Özüak et al., 2014; Sachs et al., 2015). The red flour beetle Tribolium castaneum provides an interesting intermediate position within this evolutionary series. In both Drosophila and Tribolium, Toll signalling leads to the formation of a ventral-to-dorsal nuclear gradient of the NF-κB transcription factor Dorsal. The NF-κB/Dorsal gradient is fairly stable in Drosophila and patterns the DV axis by activating or repressing target genes in a concentration-dependent manner (Hong et al., 2008; Ozdemir et al., 2014; Reeves and Stathopoulos, 2009; Reeves et al., 2012; Roth et al., 1989; Rushlow et al., 1989; Steward, 1989). Among the targets of NF-κB/Dorsal are transcription factors directly involved in cell fate specification such as twist (twi), as well as components of other signalling pathways which indirectly specify subregions along the DV axis. In particular, the NF-κB/Dorsal gradient tightly controls a number of BMP pathway components, such as the secreted BMP inhibitor short gastrulation (sog) and the BMP ligand decapentaplegic (dpp), and thereby initiates the formation of a BMP signalling gradient which specifies dorsal cell fates (Francois et al., 1994; Jazwinska et al., 1999; O'Connor et al., 2006; Wang and Ferguson, 2005). Thus, Toll signalling in Drosophila uses multiple mechanisms to provide stable patterning information for the entire DV axis.

In contrast to the situation in Drosophila, Toll signalling in Tribolium is dynamic. The NF-κB/Dorsal gradient is initiated in a broad domain that then rapidly refines and disappears while the expression domains of potential target genes are established (Chen et al., 2000; Nunes da Fonseca et al., 2008). The number of genes controlled by Toll signalling in Tribolium appears to be reduced compared with Drosophila. For example, we have not found genes that are repressed by Tc-NF-κB/Dorsal and only one mechanism, the ventral activation of Tc-sog, is employed to polarize BMP signalling (van der Zee et al., 2006). BMP signalling appears to play a larger role in DV patterning than in Drosophila (Nunes da Fonseca et al., 2010; van der Zee et al., 2006). This is even more pronounced in more basally branching insects. In the jewel wasp Nasonia vitripennis, representing the Hymenoptera (the outgroup to all other holometabolous insects), BMP patterns the entire axis, whereas Toll is only required locally to induce mesoderm formation (Özüak et al., 2014). In the milkweed bug Oncopeltus fasciatus (a hemimetabolous insect), Toll lacks all instructive roles for tissue specification; it only polarizes BMP signalling, which in turn patterns the DV axis (Sachs et al., 2015). Taken together, the DV gene regulatory network (GRN) of Tribolium deviates from that of Drosophila in many interesting ways, but has not diverged as much as the DV-GRN of more basally branching insects, making it a preferable choice for direct molecular comparisons with Drosophila.

Aside from the attractive intermediate position of its DV-GRN, Tribolium offers a number of additional advantages. First, Tribolium development shares many features with other insects (including species from basal lineages), which are absent or highly modified in Drosophila, e.g. the short-germ mode segmentation and the presence of two complete extraembryonic membranes, serosa and amnion, instead of a reduced amnioserosa (Beermann and Schröder, 2008). Second, Tribolium provides many tools for functional studies. It possesses a particularly strong systemic RNAi response – which was the basis for a recent genome-wide RNAi screen (Schmitt-Engel et al., 2015) – and techniques for forward genetic screens, transgenesis and live imaging are available (Benton et al., 2013; Berghammer et al., 1999, 2009; Koelzer et al., 2014; Sarrazin et al., 2012). Therefore, we have chosen Tribolium for unbiased genomic approaches like ChIP-seq (D. Stappert, Two novel, complementary next generation sequencing approaches to reveal the dorso-ventral gene regulatory network of Tribolium castaneum. PhD Thesis, University of Cologne, 2014) and RNA-seq, which should enable us to reconstruct the Tribolium DV-GRN for a thorough comparison with Drosophila and more basally branching insects. Here, we describe the results of RNA-seq experiments that compare transcriptomes of strongly ventralized or dorsalized embryos generated by knockdown of key DV patterning genes by parental RNAi. We also include embryos lacking only the mesoderm generated by Tc-twi RNAi and, for the first time, provide a thorough description of the Tc-twi knockdown phenotype.

We used four highly penetrant knockdown conditions to produce narrowly staged embryos with DV defects for transcriptome analysis (Fig. 1; Fig. S1). Two conditions were selected to produce severely dorsalized embryos: Tc-Toll and Tc-sog RNAi. Both lead to a loss of ventral (neuroectodermal) and an expansion of dorsal ectodermal cell fates. However, only Tc-Toll knockdown yields embryos in which the mesoderm is deleted (apparent from the loss of Tc-twi expression; Fig. 1A-C; Fig. S1H,I) (Nunes da Fonseca et al., 2008; van der Zee et al., 2006). One condition was used to produce severely ventralized embryos: Tc-dpp RNAi, which results in an expansion of the neuroectoderm at the expense of the dorsal ectoderm (Fig. 1A,D) (van der Zee et al., 2006). In Tribolium, DV and anterior-posterior (AP) patterning are more obviously linked than in Drosophila, as shown by the fact that the Tribolium head anlage occupies a ventral wedge region of the blastoderm and, similarly, the serosa is derived from an anterior-dorsal region of the blastoderm. As a result of this link, dorsalization of embryos reduces the size of the head primordium whereas ventralization causes an expansion (Fig. 1A-D). The only knockdown phenotype used in this study that does not affect AP patterning is that of Tc-twi (Fig. 1E).

Whereas a furrow forms along the ventral midline in wild-type embryos (Handel et al., 2005), Tc-twi RNAi embryos lack furrow formation (Fig. 2A,B). The tissue along the ventral midline, where Tc-twi is usually expressed, seems to be static and appears plate-like from the last mitotic division before gastrulation until germ band extension starts (Movies 1-3) (for staging, see Benton et al., 2013). The germ rudiment around this static plate condenses and the posterior amniotic fold and the head lobes form as in wild-type embryos (Fig. 2A,B; Movies 1-3). During germ band extension, Tc-twi RNAi embryos are expanded laterally and are less compact than wild-type embryos (Fig. 2C,D and Movies 1-3). The most obvious morphological difference is observed at the posterior end of fully segmented germ band embryos. In wild-type embryos, the posterior end of the germ band shows the opening of the hindgut (Fig. 2E,G) that extends into the interior and harbours the Malpighian tubules at its anterior end (King and Denholm, 2014). In Tc-twi RNAi embryos, the hindgut is detached from its usual location and flipped outwards so that it appears as an external tube connected to the posterior tip of the embryo (Fig. 2F,H). The Malpighian tubules, which are not clearly visible in unstained wild-type embryos, become finger-like protrusions at the tip of this rotated hindgut (Fig. 2E-H).

The lack of the ventral furrow through which mesoderm cells invaginate during gastrulation (Handel et al., 2005) raises the question of whether any cells invaginate and become mesoderm in the absence of Tc-twi. Cross sections through anterior regions of gastrulating wild-type embryos show the invaginating mesodermal cells that have lost epithelial organization and are clearly distinct from, and flanked by, the ectodermal epithelium (Fig. 2I). The same is observed in cross sections through posterior regions, although here the number of invaginating cells is larger than anteriorly (Fig. 2K). In sections through Tc-twi RNAi embryos, no invaginating cells that lose epithelial organization are present (Fig. 2J,L). This suggests that no mesodermal cells are specified at this stage. Interestingly, the knockdown of Tc-twi also seems to influence the structure of the posterior growth zone, as the layer of mesenchymal cells (presumptive mesoderm) that separate the germ cells from the ectoderm is absent after Tc-twi knockdown, as shown by sagittal sections through extending germ band embryos (Fig. 2M,N).

As the mesoderm seems to be missing in the Tc-twist knockdown embryos and because muscles are mesodermal derivatives, we checked whether muscles are missing after knockdown of Tc-twi. To do this, we used the muscle enhancer line pBA19 (Lorenzen et al., 2003). The knockdown of Tc-twi in this line resulted in the complete loss of body wall muscles (Fig. 2O,P). In addition, Tc-twi knockdown in a cardioblast enhancer line (Koelzer et al., 2014) reveals that the cardioblasts are also lost (Fig. S2).

Taken together, these findings suggest that upon Tc-twi knockdown, no mesodermal derivatives form from the ventral-most regions of blastoderm embryos and from the posterior growth zone during later developmental stages. Based on this data we suggest that, as in Drosophila, knockdown of twi in Tribolium leads to the absence of most, if not all mesodermal derivatives during later development. However, while Drosophila twi mutants have a twisted germ band, leading to a highly abnormal cuticle phenotype (Simpson, 1983), the loss of mesodermal derivatives has no detectable consequences for ectodermal development in Tribolium. Tc-twi RNAi embryos show largely normal germ band extension (Fig. 2H) and retraction (Fig. S2B). In addition, extraembryonic development, including dorsal closure, occurs as in the wild type, and the secreted cuticles lack overt defects (Fig. S2C,D).

We knocked down Tc-Toll, Tc-twi, Tc-dpp and Tc-sog by pRNAi and subsequently sequenced the transcriptomes of the respective knockdown embryos. We then compared the transcriptomes of knockdown and control embryos to identify genes that are regulated by Tc-Toll, Tc-twi or by BMP signalling. As a control, the transcriptomes of offspring laid by mothers injected with dsRNA specific for dsRed were sequenced. dsRed is not encoded in the Tribolium genome and should not specifically affect gene expression, but this treatment accounts for possible expression changes caused by injection of dsRNA. To investigate if there are such unspecific expression changes we also sequenced the transcriptomes of embryos derived from untreated wild-type mothers.

After verification of the general knockdown penetrance in embryos laid by injected mothers (>90% for each treatment, data not shown), RNA was prepared from narrowly staged embryos (7.5 h to 11.5 h after egg laying at 30°C), which ranged in age from undifferentiated blastoderm to early gastrulation (horseshoe amniotic fold stage; Benton et al., 2013). The knockdown efficiency (∼60-70%) in each individual RNA sample was verified by qRT-PCR (Fig. S1). Finally, 24 samples were chosen for sequencing: five biological replicates each for Tc-Toll knockdown, Tc-twi knockdown and wild-type embryos and three biological replicates each for Tc-sog, Tc-dpp and dsRed knockdown embryos. The data were processed (see Material and Methods) and tested for differentially expressed genes using two independent algorithms. Only genes that were detected to be differentially expressed by both algorithms were further considered. Interestingly, no genes were differentially expressed between untreated embryos and embryos from mothers injected with dsRed [false discovery rate (FDR) of 5%]. This is an important result in itself, as it indicates that gene expression is only affected by pRNAi against endogenous genes and shows that untreated embryos are a sufficient control in future experiments. To gain additional statistical power for analysis of the knockdown transcriptomes, transcriptomes from both dsRed pRNAi and wild-type embryos were treated as controls.

In total, 796 genes were differentially expressed with a FDR of 1% compared with control embryos (Fig. 3A): 347 upon Tc-twi, 310 upon Tc-Toll, 18 upon Tc-sog and 377 upon Tc-dpp knockdown (Tables S1-S5). To narrow down the large numbers of differentially expressed genes, we focused on those 222 genes that were differentially expressed in more than one knockdown situation (Fig. 3A; Table S6). As expected, the genes regulated by dorsalizing and ventralizing genes were predominantly regulated in opposite directions (e.g. dpp and Toll). Genes co-regulated by two ventralizing genes were mostly regulated in the same direction (e.g. Toll and twist) (Fig. 1 and Fig. 3B).

To find genes that act in the neuroectoderm or dorsal ectoderm, we focused on comparison of Tc-dpp and Tc-sog knockdown. As described earlier, Tc-dpp knockdown leads to a massive expansion of neuroectoderm, while Tc-sog knockdown results in the opposite phenotype, that being an almost complete loss of the neuroectoderm. Both knockdowns have little impact on the mesoderm (Fig. 1). We found 14 genes to be differentially expressed in these opposing situations at an FDR of 1% (Fig. 3A; Table S7). Of these, 10 genes are downregulated upon Tc-dpp knockdown and upregulated upon Tc-sog knockdown (Fig. 3B; Table S7) and the remaining 4 genes show the reverse behaviour (Fig. 3B). As predicted from the fate map shifts (Fig. 1), genes belonging to the former group are mainly, and often exclusively, expressed in the serosa (Fig. 4A-E, data not shown for TC006222, TC010855, TC001439 and TC015188; TC008855 could not be cloned), while genes belonging to the latter group are expressed in the neurogenic ectoderm (Fig. 4F-H). The neuronal cell fate determinant TC010596/Tc-prospero (Fig. 4H and data not shown) (Biffar and Stollewerk, 2014; Ungerer et al., 2011) and the leucine-rich repeat (LRR) transmembrane protein TC014658/Tc-tartan (Fig. 4G and Fig. 5F) (Chang et al., 1993) both show broad neuroectodermal expression along the entire AP axis. In addition, TC014658/Tc-tartan is upregulated in a wedge-shaped region of the head anlage. The other two genes, TC004745/Tc-patched (Fig. 4F) (Farzana and Brown, 2008) and TC007409/Tc-twin of eyeless (data not shown, see Yang et al., 2009), show high expression in the head region of differentiated blastoderm embryos. This region is lost upon dorsalization and expanded upon ventralization (Fig. 1). To search the dataset for more candidates we also included some genes selected with a FDR of 5% (Table S7), e.g. TC000871/Tc-uninflatable (Tc-uif), which is known to be involved in Notch signalling in Drosophila (Xie et al., 2012; Zhang and Ward, 2009).

We next functionally tested three of the identified genes: Tc-uif, Tc-tartan and one of the genes expressed in the serosa, TC006771. For each knockdown, we analysed survival of the injected mothers, the morphology of DAPI-stained embryos and the expression of marker genes for DV cell specification: Tc-pannier (Tc-pnr, amnion; van der Zee et al., 2005), Tc-twi (mesoderm) and Tc-achaete-scute homolog (Tc-ash, neuronal precursors; Wheeler et al., 2003). Tc-engrailed (Tc-en) and/or Tc-gooseberry (Tc-gsb) (Brown et al., 1994; Davis et al., 2001; Nunes da Fonseca et al., 2008) (Fig. S3) were used to monitor segmentation.

Tc-uif is first expressed in the neurogenic ectoderm with increasing expression levels towards the border of the presumptive mesoderm (Fig. 5A,B). Later, Tc-uif is expressed in the whole germ band with the exception of the mesoderm (Fig. 5C). Tc-uif knockdown leads to embryos with thinner and twisted abdomens (61.45%, n=182) (Fig. 5D,E). However, none of the marker genes used to analyse the phenotype showed altered expression patterns in early embryos (Fig. S3). Thus, despite the early expression of Tc-uif, the phenotype apparently arises during later development within the growth zone.

Early neuroectodermal and head expression of Tc-tartan (Fig. 4G and Fig. 5F) is followed by segmental expression in later stages (Fig. 5G,H). Tc-tartan is not expressed in the mesoderm; however, upon knockdown of Tc-twi (Fig. 5I, compare with Fig. 4G), its expression expands into the future mesoderm, suggesting a suppressive function of Tc-twi or of a Tc-Twi target gene. Expression analysis of Tc-ash shows that Tc-tartan is involved in neurogenesis because cells expressing Tc-ash (marking proneural clusters) were largely lost upon knockdown of Tc-tartan (Fig. 5J,K).

Upon TC006771 knockdown, we observed embryos with interrupted germ bands and thin abdominal segments. However, the penetrance of this phenotype was low (6.86% of all embryos, n=322). Although comparable defects are not seen in controls, we did not further investigate this phenotype. Taken together, by comparing the transcriptomes of severely ventralized embryos with that of severely dorsalized embryos, we identified new genes expressed in the serosa and discovered new early functions for genes that are only required at later stages in Drosophila (see Discussion). In particular, Tc-tartan appears to be a component of the Tribolium DV-GRN required for early neurogenesis.

While Tc-Toll knockdown results in a strong dorsalization, Tc-twi knockdown has only minor effects on the early fate map (Fig. 1). Furthermore, Tc-twi is downstream of Tc-Toll (Nunes da Fonseca et al., 2008). Therefore, we hoped to identify more direct targets of the early DV-GRN by focusing on the group of genes that are downregulated upon knockdown of both Tc-Toll and Tc-twi. This assumption is supported by the observation that both Tc-twi and Tc-cactus (which is dependent on Tc-Toll and Tc-twi) are among the 19 genes that fall into this category with an FDR of 1%. (Nunes da Fonseca et al., 2008) (Table S8; Fig. 6A,B). Further, Tc-snail seems to be regulated by Tc-Toll and Tc-twi (Fig. 6C), like in Drosophila (Reeves and Stathopoulos, 2009). The same applies for the homologues of the Drosophila twi target genes heartless (htl) and down of FGF receptor (dof) (Beermann and Schröder, 2008; Beiman et al., 1996; Gisselbrecht et al., 1996; Imam et al., 1999; Vincent et al., 1998) (Fig. 6D,E). Interestingly, Tc-Delta expression depends on Tc-Toll and Tc-twi and is found in the presumptive mesoderm in wild-type embryos (Fig. 6F). This is strikingly different from Drosophila, where Delta is repressed in the presumptive mesoderm and shows graded expression in the neurogenic ectoderm (Vassin et al., 1987). Along with Tc-Delta, the two enhancer of split homologues of Tribolium Tc-E(spl)1 and Tc-E(spl)3 (Kux et al., 2013) are regulated by Tc-Toll and Tc-twi and are expressed in a mesodermal domain with higher levels toward the lateral borders (Fig. 6G,H). Expression of E(spl) complex genes within the presumptive mesoderm has also not been reported for Drosophila (Knust et al., 1992; Wech et al., 1999).

The genes TC001667/integrin αPS2 (inflated), TC005328, TC010105/l(2)efl and TC010461 were also downregulated upon Tc-twi and Tc-Toll knockdown, and are weakly expressed in the mesoderm of wild-type embryos (Fig. 6I-L). The remaining genes from this group (TC003461, TC003606, TC007056, TC009862, TC010195, TC013142) were cloned, but did not show a specific expression pattern, with the exception of the pair-rule gene TC008064/Tc-sloppy paired (Choe and Brown, 2007), which is expressed in stripes along the anterior-posterior axis (Fig. 6M).

In summary, most of the genes that are downregulated upon knockdown of Tc-twi and Tc-Toll are expressed in the presumptive mesoderm shortly before gastrulation. Given their early expression pattern they might be directly regulated by both Tc-Toll and Tc-twist, or by Tc-twist alone.

For functional analysis, we chose TC004713/Tc-htl, the sole FGF receptor of Tribolium and TC011323/Tc-dof, a cytoplasmic signal transducer downstream of the FGF receptor (Fig. S4), as well as two genes selected with reduced stringency (5% FDR): TC011114/Tc-zinc-finger homeodomain 1 (Tc-zfh1) and TC005184/Tc-LamininB1 (Tc-LanB1, Fig. S4). All four genes are expressed in the presumptive mesoderm at the primitive pit stage (Fig. 6D,E and Fig. 7E,K).

Knockdown of Tc-htl or Tc-dof does not affect expression of the aforementioned marker genes suggesting no major patterning defects during germ band extension (Fig. S4A-M). However, using the transgenic line G04609 (Koelzer et al., 2014), we observed that knockdown of Tc-htl resulted in complete loss of cardioblasts (Fig. 7A,B) confirming earlier observations (Sharma et al., 2015), whereas Tc-dof knockdown resulted in a partial loss of this tissue (Fig. 7A,C). In addition, Tc-htl embryos show a mispositioned hindgut with Malpighian tubules protruding at the tip, as in Tc-twi knockdown embryos (compare Fig. 7B with Fig. 2F and Fig. S2B,D; see also Sharma et al., 2015). This hindgut phenotype was, however, never observed in Tc-dof knockdown embryos (Fig. 7C,D). In the course of a genome-wide RNAi screen in Tribolium, the impact on somatic muscle differentiation was also significantly weaker for knockdown of Tc-dof than Tc-htl (Donitz et al., 2015; Schmitt-Engel et al., 2015). This indicates that Tc-dof is not required for all aspects of FGF signalling in Tribolium, which is in contrast to Drosophila where dof mutants almost completely abolish FGF signalling (Imam et al., 1999; Vincent et al., 1998). The fact that Tc-htl and Tc-twi knockdown cause the same hindgut phenotype suggests that Tc-twi primarily acts via FGF signalling to control hindgut morphogenesis.

Tc-zfh1 expression is maintained in the mesoderm during gastrulation and germ band extension (Fig. 7F,G). Knockdown of Tc-zfh1 results in reduced survival of injected adult beetles (28% survival, n=200; compared with Tc-Toll and Tc-twist injected beetles with 72% and 78% survival, respectively, n=300) and reduced egg production. Reduction of the dsRNA from 1 µg/µl to 0.1 µg/µl alleviated the sterility effect and allowed expression analysis of marker genes (Fig. S4N-Q). Interestingly, of the four investigated marker genes, only Tc-twi expression was affected. Although early Tc-twi expression is normal upon knockdown of Tc-zfh1, Tc-twi is lost in abdominal segments during germ band extension (Fig. 7, compare I with J). This suggests that Tc-zfh1 is required to maintain Tc-twi expression subsequent to the activation of Tc-zfh1 by Tc-twi. In Drosophila, mutations of zfh1 only cause weak muscle phenotypes (see Discussion) (Lai et al., 1993).

Tc-LanB1 is expressed in the entire germ rudiment and upregulated in the presumptive mesoderm (Fig. 6N and Fig. 7K). Tc-LanB1 knockdown embryos show no early embryonic defects based on the selected marker genes (Fig. S4R-U). However, the abdomen of extended germ band embryos is thinner and irregularly shaped. Tc-LanB1 knockdown in the muscle enhancer line pBA19 reveals a broad spectrum of muscle defects ranging from mild disorganization to complete loss of body wall muscles (Fig. 7M and Fig. S5A-D), whereas cardioblasts are always lost (Fig. S5E-H). In Drosophila, the complete loss of LanB1 function leads to much milder defects during late organogenesis (Urbano et al., 2009).

Taken together, selecting genes that were downregulated after Tc-Toll and Tc-twi RNAi led to the discovery of new expression patterns, such as the mesodermal expression of the Tc-E(spl) genes, as well as new early phenotypes caused by the knockdown of genes known to have later roles in Drosophila. In particular, Tc-zfh1 is likely to be a component of the Tribolium DV-GRN that is required to maintain the expression of Tc-twi, the key regulator of mesoderm development.

To understand the Tribolium DV-GRN it is necessary to identify all of its components, irrespective of homologues in other species. By knocking down key regulators of DV patterning followed by comprehensive transcriptome analysis, we identified over 750 differentially expressed genes. This dataset is a valuable resource, not only for the identification of early, locally expressed genes, but also for the unbiased identification of DV network components. In this study, we focused only on those genes that were differentially expressed in more than one knockdown situation. Even this resulted in a number (222) that precluded complete and careful analysis. We therefore selected two subgroups for in-depth study: the 14 genes differentially expressed between Tc-dpp, Tc-sog and Tc-Toll knockdowns and the 19 genes downregulated upon Tc-Toll and Tc-twi knockdown. We expected to identify ectodermal genes in the first group and mesodermal genes in the second group. Five genes from the first group show expression only in the (anterior-dorsal) serosa and three are upregulated in the (anterior-ventral) head region. This high number of genes with differential expression along the AP axis reflects the aforementioned fact that the sizes of the serosa and head anlagen are affected by global DV patterning genes like Tc-Toll, Tc-sog and Tc-dpp. In contrast, the second group requiring downregulation following Tc-twi knockdown did not result in the selection of AP patterning genes, and thus is more likely to be enriched in bona fide elements of the DV-GRN. However, further experiments such as ChIP-seq or detailed molecular studies will be necessary to understand the regulatory wiring connecting these genes. From our data, we chose a limited number of genes for functional analysis by RNAi. In the following, we will first discuss the Tc-twi knockdown phenotype before returning to the genes identified through comparative transcriptome analysis.

Live imaging of knockdown embryos in enhancer trap lines demonstrated that Tc-twi function is largely confined to the mesoderm and that the loss of the mesoderm has little effect on the rest of the embryo. With the exception of ventral furrow formation, all early (folding of the amnion, germ band extension) and most of the late morphogenetic movements (germ band retraction, opening of the amniotic cavity and dorsal closure) are not obviously affected (Fig. 2; Fig. S2A,B). Only the hindgut is not positioned correctly. The Tc-twi knockdown embryos secrete a cuticle that is hard to distinguish from controls (Fig. S2C,D). However, Tc-twi knockdown larvae completely lack somatic mesoderm as well as cardioblasts. Thus, in Tribolium, ectodermal development occurs apparently largely independent of mesodermal specification. This phenotype is remarkable as the Drosophila twi mutants have been identified by virtue of their strong cuticular defects, which include head defects, incomplete germ band retraction and an overall twisted appearance (Simpson, 1983).

The Tc-twi knockdown phenotype also provides interesting insights into mesoderm formation in the growth zone. In an earlier study, it was suggested that the inner mesenchymal cell layer of the growth zone is derived from cells that had been internalized through ventral furrow formation in the early embryo, which maintain a mesodermal identity and later give rise to the mesoderm of the newly forming segments (Handel et al., 2005). However, the lack of Tc-twi expression in these inner layer cells represented a caveat for this hypothesis, because in Drosophila, mesodermal precursor cells continue to express twi. The Tc-twi knockdown phenotype resolves this problem, as the growth zone of the knockdown embryos lack the inner mesenchymal cell layer. Thus, these cells indeed are derived from cells internalized through prior ventral furrow formation – the process that is blocked upon Tc-twi knockdown.

The hindgut phenotype of Tc-twi knockdown represents an interesting deviation from Drosophila. After Tc-twi knockdown, the hindgut does not point to the interior of the embryo where it is usually connected to the midgut, but it rather points posteriorly and the Malphigian tubules form finger-like protrusions at its tip (Fig. 2E-H). We observed a virtually identical phenotype after knockdown of Tc-htl (Fig. 7B). Tc-htl knockdown embryos lack Tc-twi-expressing cells surrounding the hindgut and are likely to give rise to the caudal visceral mesoderm (CVM) (King and Denholm, 2014; Sharma et al., 2015). Thus, in Tribolium, Tc-twi appears to be required for the formation of the CVM and this effect is probably not only indirect through regulation of FGF signalling by Tc-twi, but also direct since the CVM precursor cells express Tc-twi. In contrast, in Drosophila, the CVM cells do not express twi and represent the only group of mesodermal cells that can form in the absence of twi (Kusch and Reuter, 1999). This observation indicates that the extent to which twi was recruited to specify different mesodermal subtypes varies between insect lineages.

Among the ectodermally expressed genes, Tc-uif and Tc-tartan showed interesting early expression patterns and were selected for functional analysis. In Drosophila, uif is required during late embryogenesis (tracheal inflation; Zhang and Ward, 2009) and during imaginal disc development where it acts as a Notch antagonist (Loubéry et al., 2014; Xie et al., 2012). No early developmental functions of uif have been discovered, although the expression of Drosophila uif already starts before gastrulation in a pattern reminiscent to that of Tc-uif (Zhang and Ward, 2009). In Tribolium, uif is clearly required for embryonic development because Tc-uif knockdown affects growth zone patterning or morphogenesis.

Tartan and its close homologue Capricious are transmembrane proteins with extracellular leucine-rich repeats (LRRs). In Drosophila, tartan is first expressed in a double-segmental (pair-rule like) pattern and later in proneural clusters and sensory mother cells (Chang et al., 1993). For capricious, only late embryonic expression in motorneurons and muscles was reported. In Tribolium, we identified only one gene with similarity to both tartan and capricious. Since this Tribolium gene shares some aspects of early embryonic expression with Drosophila tartan, we named it Tc-tartan (Fig. 5G-H). Loss of tartan in Drosophila leads to weak defects during late peripheral and central nervous system development (Chang et al., 1993). As there is no information on the embryonic phenotype of the double mutant between tartan and capricious we cannot exclude that deletion of both genes has a more severe embryonic phenotype in Drosophila. In wing imaginal discs, both genes act together to stabilize the DV compartment boundary (Milán et al., 2005, 2001). However, the lack of a prominent early embryonic expression of capricious does not support a similar interaction for the early Drosophila embryo. In Tribolium, loss of Tc-tartan leads to a severe and interesting phenotype: the absence of proneural clusters marked by Tc-ash expression (Fig. 5J). This indicates that tartan has an essential function in early CNS development in Tribolium. Thus, as in the case of Tc-uif, the Tc-tartan phenotype appears to be substantially stronger than that of the corresponding Drosophila gene.

As Tc-twist knockdown only affects the mesoderm, the expression of genes that are downregulated upon both Tc-Toll and Tc-twist knockdown should be restricted to the mesoderm. This prediction was confirmed: from 19 genes selected with an FDR of 1%, 12 showed mesodermal expression. By relaxing the FDR to 5%, additional ventrally expressed genes were identified. Besides known candidates (Tc-twi, Tc-snail, Tc-cactus), representatives of two signalling pathways (Notch and FGF) and components involved in cell adhesion (integrins, laminin and dystroglycan) were identified (Table S8). Only one of the newly identified genes, Tc-zfh-1, encodes a transcription factor.

While expression and function of the FGF pathway components Tc-htl (the single FGF receptor) and Tc-dof (an adaptor protein) is very similar to that of the Drosophila homologues (Beermann and Schröder, 2008; Sharma et al., 2015), the expression of Notch pathway components deviates strikingly. In particular the broad ventral expression of the Notch ligand Delta (Aranda et al., 2008) and the two Tribolium E(spl) homologues suggest that Notch signalling is active in the presumptive mesoderm. In Drosophila, Notch signalling is inhibited in the mesoderm through several parallel mechanisms, including the repression of Delta (Bardin and Schweisguth, 2006; De Renzis et al., 2006; Vassin et al., 1987).

The early mesodermal expression of the extracellular matrix receptor αPS2 integrin (inflated) is shared between Drosophila and Tribolium (Fig. 6I) (Bogaert et al., 1987; Wehrli et al., 1993). In Tribolium, but not Drosophila, however, the genes encoding the basement-membrane component Laminin B1 and Dystroglycan (another ECM receptor) are upregulated in the presumptive mesoderm (Fig. 6N and Fig. 7K; Table S8, data not shown) (Dekkers et al., 2004; Montell and Goodman, 1989; Schneider and Baumgartner, 2008; Urbano et al., 2009). The interaction of laminin with surface receptors like integrins and dystroglycan plays an important role in gastrulation in many organisms, in particular in vertebrates (Ettensohn and Winkelbauer, 2004). Thus, the early expression of these components in Tribolium might indicate an ancestral feature of ventral furrow formation which has been partially lost in Drosophila. Although the Tc-LanB1 knockdown does not abolish ventral furrow formation, it might influence its kinetics. The earliest phenotypes we observed are morphological defects in the posterior germ band. In late embryos, the body wall muscles are highly disorganized (Fig. S5). Thus, in Tribolium, LanB1 is required earlier than in Drosophila for correct morphogenesis and the late embryonic knockdown phenotypes appear to be much stronger than the LanB1 mutant phenotypes of Drosophila, which show rather mild muscle defects (Urbano et al., 2009).

One of the most interesting genes we have found is Tc-zfh1. This gene had also been identified through differential transcriptome analysis in Drosophila (Casal and Leptin, 1996). It was first described because it codes for an unusual transcription factor combining zinc fingers and a homeodomain (Casal and Leptin, 1996; Fortini et al., 1991; Lai et al., 1991). The Drosophila zfh1 mutants show weak somatic muscle phenotypes (Lai et al., 1993), lack the CVM and have defects in the gonadal mesoderm (Broihier et al., 1998; Kusch and Reuter, 1999). Tc-zfh1 knockdown leads to multiple morphogenetic defects during germ band formation. Most importantly, however, Tc-zfh1 is required for maintaining Tc-twi expression and thus is an essential component for mesoderm development in Tribolium.

In conclusion, by studying subsets of genes differentially expressed after manipulating DV patterning in Tribolium we have succeeded in identifying components of the Tribolium DV-GRN in an unbiased way. The dataset presented here will be screened for further candidates and will serve as a strong foundation for more detailed studies on a wide range of questions surrounding insect embryo axis patterning, cell fate specification and morphogenesis.

The San Bernadino wild-type strain was used for RNAi injections and as wild-type control, if not specified otherwise. For selected experiments, we used pBA19 (Lorenzen et al., 2003), G04609 (Koelzer et al., 2014) and an nGFP strain (Sarrazin et al., 2012).

ISH was essentially performed as described previously (van der Zee et al., 2005) and the complete protocol can be found in supplementary Materials and Methods.

dsRNA was produced with the MEGAscript T7 Kit (Ambion) according to the manufacturer's protocol. The synthesized dsRNA was purified by phenol chloroform extraction. For all knockdown experiments, with the exception of the Tc-dpp knockdown, female pupae were injected with 1 µg/µl dsRNA. For the Tc-dpp knockdown adult females were injected. For TC011114, phenotypes were too strong after injection of 1 µg/µl and we injected 0.5 µg/µl and 0.1 µg/µl instead. For off-target analysis of primers and RNA fragments used to produce dsRNA, see Table S10.

Total RNA was isolated using TRIzol reagent (Invitrogen) from a batch of embryos staged to 7.5 h to 11.5 h after egg lay at 30°C. See supplementary Materials and Methods for a step-by-step protocol of RNA isolation. The knockdown efficiency was evaluated by qRT-PCR for each individual RNA sample (e.g. see Fig. S1). Five biological replicates for each Tc-Toll knockdown, Tc-twist knockdown and wild-type embryos and three biological replicates for each Tc-sog, Tc-dpp and dsRed knockdown embryos were sequenced in 100 bp mode on a HiSeq 2000 (Illumina) by the Cologne Center for Genomics (CCG, paired-end mode) or by GATC Biotech (Konstanz, single-end mode). Biases introduced by sequencing different samples in different sequencing facilities were detected by a principle component analysis (data not shown) and considered in all downstream analyses. Quality of the FASTQ formatted sequences was verified using FastQC (v.3.5.12; www.bioinformatics.babraham.ac.uk/projects/fastqc/). Bowtie2 (v.2.0.2; bowtie-bio.sourceforge.net/bowtie2) with ‘very sensitive’ settings used in Tophat (ccb.jhu.edu/software/tophat/) to map the sequences to the Tcas3.0 gene set. The data have been submitted to NCBI under accession number PRJNA315762. Count tables that summarize reads per gene were compiled in R (v.3.0.2 2013-09-25, accessed with RStudio v.0.97.551) using either the ‘countOverlaps’ function from the IRanges (v. ≥1.18.2; bioconductor.org/packages/release/bioc/html/IRanges.html) package (for subsequent processing with edgeR) or the ‘summarizedOverlaps’ function in mode ‘union’ from GenomicRanges (v.1.12.4; bioconductor.org/packages/release/bioc/html/GenomicRanges.html) (for subsequent processing with DESeq2). The data were then tested for differentially expressed genes with edgeR (v.3.2.4; bioconductor.org/packages/release/bioc/html/edgeR.html) and DESeq2 (v.1.0.18; bioconductor.org/packages/release/bioc/html/DESeq2.html) (Anders and Huber, 2010; Robinson et al., 2010). The edgeR analysis was essentially performed as described in section 4.5 of the edgeR Users Guide (31 March 2013). The DESeq analysis was performed as described in the Analyzing RNA-Seq data with the DESeq2 package manual (15 July 2013), with the exception that the function ‘nbinomWaldTest’ was used without using a cut-off for Cook's distance. Only genes that were detected to be differentially expressed by both algorithms at a 5% FDR were considered to be differentially expressed.

The efficiency of RNAi treatment was evaluated with qRT-PCR prior to transcriptome sequencing. A part of each staged RNA batch isolated from knockdown and wild-type embryos was transcribed into cDNA using SuperScriptVILO cDNA Synthesis Kit (Life Technologies). Quantities of the knocked down genes or their target genes were compared between offspring from dsRNA-injected mothers and wild-type embryos using SYBR Green (Applied Biosystems) and a 7500 Fast Real-time PCR system (Applied Biosystems). The reaction volume was 25 µl (12.5 μl SYBR Green PCR master mix, 2 μl cDNA, 1 μl forward primer, 1 μl reverse primer, 8.5 μl H₂O) and the cycling profile was: 94°C 5 min, 35 cycles of 94°C 30 s, 55°C 30 s, 72°C 1 min and extension at 72°C for 5 min. The amplification efficiency for each primer pair was empirically determined in each qRT-PCR run and for calculation of the RNA ratios from knockdown embryos and wild-type embryos the ΔCt method was used. See primer list for primer sequences.

Live imaging was performed with dechorionated embryos mounted on slides and immersed in Halocarbon oil 700 (Sigma) on an Applied Precision DeltaVision RT wide-field microscope at 30°C. Movies were generated using an ImageJ macro created by Thorsten Horn (Institute of Developmental Biology, University of Cologne).

DAPI and Phalloidin staining were carried out using standard methods as detailed in supplementary Materials and Methods. Antibody staining for phosphotyrosine was performed using 1:500 monoclonal anti-phosphotyrosine (Sigma, P3300) and 1:400 goat anti-mouse secondary antibody, Alexa Fluor 555 conjugate (Thermo Fisher Scientific, A-21422).

Full list of primers used is shown in Table S9.

We are especially grateful to Rodrigo Nunes da Fonseca and Jeremy A. Lynch for advice and discussions. We thank Thorsten Horn for his help with live imaging, and Matt Benton, Matthias Pechmann and Waldemar Wojciech for valuable comments on the manuscript and for corrections.