|Home | About | Journals | Submit | Contact Us | Français|
After fertilization the embryonic genome is inactive until transcription is initiated during the maternal-zygotic transition1,2,3. This transition coincides with the formation of pluripotent cells, which in mammals can be used to generate embryonic stem cells. To study the changes in chromatin structure that accompany pluripotency and genome activation, we mapped the genomic locations of histone H3 molecules bearing Lysine trimethylation modifications before and after the maternal-zygotic transition in zebrafish. Trimethylation of Lysine 27, which is repressive, and trimethylation of Lysine 4, which is activating, were not detected before the transition. After genome activation, more than 80% of genes were marked by Lysine 4 trimethylation, including many inactive developmental regulatory genes that were also marked by Lysine 27 trimethylation. Sequential chromatin immunoprecipitation demonstrated that the same promoter regions had both trimethylation marks. Such bivalent chromatin domains also exist in embryonic stem cells and are thought to poise genes for activation while keeping them repressed4,5,6,7,8. In addition, we found many inactive genes that were uniquely marked by Lysine 4 trimethylation. Despite this activating modification, these monovalent genes were neither expressed nor stably bound by RNA polymerase II. Inspection of published datasets revealed similar monovalent domains in embryonic stem cells. Moreover, Lysine 4 trimethylation marks could form in the absence of both sequence-specific transcriptional activators and stable association of RNA pol II, as indicated by the analysis of an inducible transgene. These results suggest that bivalent and monovalent domains might poise embryonic genes for activation and that the chromatin profile associated with pluripotency is established during the maternal-zygotic transition.
After fertilization, animals go through cleavage divisions that transform the 1-cell egg into a multicellular embryo. During this phase, the genome is often inactive, and zygotic transcription is only initiated during the maternal-zygotic transition. In zebrafish, the process of zygotic transcription begins three hours post fertilization (hpf), when a period of rapid synchronous cleavages ends and cell cycles lengthen9. At this time, embryonic cells (blastomeres) are pluripotent and largely identical9. To study the changes in chromatin structure that accompany pluripotency and zygotic genome activation, we mapped the genomic locations of histone H3 modifications and RNA polymerase II (RNA pol II) during this transition. Previous studies have established that nucleosomes with histone H3 Lysine 4 trimethylation (H3K4me3, catalyzed by Trithorax group proteins) and histone H3 Lysine 27 trimethylation (H3K27me3, catalyzed by Polycomb Group proteins) are associated with activation and repression of transcription, respectively10. Nucleosomes with histone H3 Lysine 36 trimethylation (H3K36me3) mark regions of transcriptional elongation11. By performing Western blot analysis on embryos from early blastula stages to gastrulation, we found that RNA pol II phosphorylation and H3K4me3 appeared during the onset of zygotic transcription and were followed by the emergence of H3K36me3 (Figure 1a). To determine the genomic locations of RNA pol II, H3K36me3, H3K4me3 and H3K27me3 before and after the maternal-zygotic transition, we performed chromatin immunoprecipitation (ChIP) and identified coprecipitated DNA fragments by hybridization to custom-designed NimbleGen tiling microarrays12 covering ~31 Mb of the zebrafish genome. Peaks were detected using MA2C13, a method that normalizes log2 fold enrichment by probe sequence and signal distribution between ChIP and control, and assigns each probe a normalized, window-smoothed MA2Cscore. Enriched domains were then defined and associated with RefSeq genes as described in the Methods section and Supplementary Figure 1. Consistent with the distribution profiles found in cultured cells14, RNA pol II was enriched near transcription start sites (Figure 1b), and H3K36me3 was detected in gene bodies, peaking towards the 3′ end (Figure 1c). H3K4me3 was localized to promoter regions (Figure 1d), while H3K27me3 covered promoters or entire genes (Figure 1e, Supplementary Figure 2). Strikingly, the comparison of ChIP-chip profiles prior to and after the maternal-zygotic transition (Figure 1b–e and Supplementary Figure 3) revealed that RNA pol II and the analyzed histone modifications were not associated with the genome before zygotic genome activation. Inspection of specific genes and genomic regions confirmed the lack of H3 methylation marks before the maternal-zygotic transition as well as their presence thereafter (Figure 1f and Supplementary Figures 4 and 5). These results indicate that RNA pol II and the histone modifications H3K36me3, H3K4me3 and H3K27me3 become associated with genomic loci only during the onset of zygotic genome activation.
We next analyzed the histone modification patterns that emerged after the maternal-zygotic transition (Supplementary Figure 6 and Supplementary Table 1; replicate analysis in Supplementary Figure 7). Of 683 RefSeq genes, 82% were occupied by H3K4me3 at their promoter. These included 97% of the 367 genes that were zygotically expressed and 64% of the 316 genes that were inactive (i.e. genes not occupied by H3K36me3, see Supplementary Discussion 1 and Supplementary Figures 1 and 6). The inactive genes were of particular interest, because studies in ES cells found that developmental regulatory genes are often inactive and associated with H3K4me3/H3K27me3 domains4,5,6,7,8. We found that in zebrafish blastomeres 36% of inactive genes were associated with H3K4me3 and H3K27me3 (Figure 2a, b and Supplementary Figure 6). Gene Ontology (GO) analysis showed that these inactive genes were enriched for those annotated as ‘DNA binding’, ‘transcription factor activity’ and ‘multicellular organismal development’, with p-values of 3.3×10−13, 8.9×10−12 and 4.1×10−6, respectively, similar to that seen for ES cells4,6,7,8. Among the H3K4me3/H3K27me3 genes were the Hox clusters (Supplementary Figure 4), cmyb (Figure 2b), evx1, and irx3. These results indicate that H3K4me3 and H3K27me3 mark regulatory genes in both mammalian ES cells and zebrafish blastomeres.
The presence of H3K4me3 and H3K27me3 on a subset of genes might indicate the simultaneous association of both chromatin marks with a given promoter in the same cells. Such “bivalent” domains were found in ES cells and have been postulated to poise genes for activation (H3K4me3) while keeping them repressed (H3K27me3)4. Alternatively, H3K4me3 marks might occupy a given promoter in only a subset of cells, whereas H3K27me3 marks might be present in a different subpopulation. Indeed, recent experiments in Xenopus have suggested that H3K4me3 and H3K27me3 do not co-occupy the same promoters15. In that study, H3K4me3 marks were found at blastula stages but, in contrast to our results (Figure 2a, b and Supplementary Figures 4, 5 and 6), H3K27me3 marks were only detected after the onset of gastrulation. This temporal separation of histone modifications would preclude bivalency in Xenopus blastomeres. This study also identified few, if any, bivalently marked promoters during gastrulation15. To determine if bivalent domains exist in zebrafish blastomeres, we performed sequential ChIP analysis (Supplementary Figure 8 and Supplementary Discussion 2). Using this assay, we found that H3K4me3/H3K27me3 genes such as cbfb, runx3 and six3b are truly bivalent (Figure 2c). These results establish that bivalent chromatin domains are present in pluripotent embryonic cells in vivo.
In addition to the 36% bivalent genes, 28% of inactive genes were associated with H3K4me3 but not H3K27me3 marks (Figure 2a, Figure 3a, b). Examples of such monovalent genes (H3K4me3 only) included psat1 (phosphoserine aminotransferase 1), which is ultimately expressed in the central nervous system and the eye, and zgc:110784, which is expressed later in the neural crest (Figure 3a). GO term analysis did not reveal any enrichment within this class. The existence of such monovalent inactive genes was surprising because H3K27me3 had been postulated to prevent transcription from H3K4me3-bound promoters4,5,16. However, closer analysis of available human and mouse ES cell data showed that many non-expressed genes also carried monovalent H3K4me3 marks (Supplementary Discussion 3 and Supplementary Figure 9). In contrast to the large fraction of monovalent H3K4me3 genes, very few genes carried monovalent H3K27me3 marks in zebrafish blastomeres (Supplementary Figure 6) and ES cells6. These results reveal remarkably similar chromatin profiles between embryonic blastomeres and embryonic stem cells.
Recent studies have indicated that genes can be poised for expression by association of RNA pol II with transcription start sites17,18,19,20,21. By examining RNA pol II occupancy, we found that only 17% of monovalent transcription units had detectable RNA pol II associated with their start sites (Figure 3b, Supplementary Figure 6). Comparing the ChIP signals for H3K4me3 and RNA pol II (Figure 3c) showed that even monovalent genes with high levels of H3K4me3 had undetectable levels of RNA pol II. Similarly, many genes in ES cells are associated with H3K4me3 but not RNA pol II (Supplementary Discussion 4). Taken together, these results reveal that many genes are initially inactive after genome activation yet marked with H3K4me3.
The presence of monovalent genes (H3K4me3 only) without detectable RNA pol II raised the possibility that sequence-specific transcription factors might direct the methylation of H3K4 before recruitment of RNA pol II22,23. Alternatively, H3K4me3 might mark monovalent genes even before the binding of sequence-specific transcription factors and recruitment of RNA pol II. Because only a subset of transcription factors is present at the maternal-zygotic transition, the latter possibility is more likely. To distinguish between these models, we analyzed an inducible transgene that encodes EGFP under the control of the heterologous transcription factor GAL4-VP16. The transgene contained 14 upstream binding sites for GAL4-VP16 joined to a minimal viral promoter (Figure 4a). Expression of the transgene was dependent on GAL4-VP16 and was not detected before zygotic genome activation (data not shown). This experimental system allowed us to compare occupancy by RNA pol II and H3K4me3 in the presence or absence of GAL4-VP16. As expected, injection of GAL4-VP16 resulted in the occupancy of the reporter gene by both H3K4me3 and RNA pol II after the maternal-zygotic transition (Figures 4b, c). In the absence of GAL4-VP16, H3K4me3 still marked the promoter, but neither RNA pol II nor H3K27me3 were detected (Figure 4b, c and data not shown). These results suggest that monovalent H3K4me3 marks can be established in the absence of sequence-specific activators and without the stable association of RNA pol II.
The chromatin profiles before and after the maternal-zygotic transition lead to three major conclusions about genome activation and pluripotency. First, genome activation is accompanied by the appearance of specific H3 trimethylation marks, suggesting that these histone modifications begin to provide regulatory information only during the maternal-zygotic transition. This finding is pertinent to the suggestion that histone marks associated with developmental regulators in human sperm prepare genes for embryonic activation or repression24. Our results suggest that if such modifications exist in zebrafish sperm, they are replaced by other marks, or they are erased or diluted during the rapid cleavage divisions before the maternal-zygotic transition. Our data are consistent with a report that identified H3K4me3 marks during the maternal-zygotic transition in Xenopus embryos15, but in contrast to that study, we detect H3K27me3 modifications already in blastomeres.
Second, our study provides evidence in embryos for the bivalent domains found in embryonic stem cells4,5,6,7,8. This finding reveals a shared chromatin landscape between transiently pluripotent cells in the embryo and permanently pluripotent cells in culture4,5. In contrast to the recent proposal that bivalent domains might not exist in embryos15, we observe clear evidence for bivalency. This discrepancy could be due to species-specific modes of gene regulation but might also be caused by experimental differences (Supplementary Discussions 2 and 5). The function of bivalent domains has not been fully established. In the context of ES cells, such marks might poise genes for activation, maintain gene repression and pluripotency and control lineage commitment25. In the developing embryo, where pluripotency is short-lived, such marks might also ensure the proper timing and quantity of developmental regulators26.
Third, we find that many genes are monovalently marked by H3K4me3 but inactive. This subset of genes has received little attention in ES cell studies, but our analysis of published datasets suggests that many non-expressed genes in ES cells also carry monovalent H3K4me3 marks. Most of these H3K4me3 domains are not associated with detectable levels of RNA pol II, and our experiments with an inducible transgene suggest that H3K4me3 occupancy can be independent of sequence-specific transcriptional activators or stable association with RNA pol II. These H3K4me3 domains might poise genes for activation by creating a platform for the transcriptional machinery27. Consistent with this idea, basal transcription factor TFIID has a higher affinity for H3K4me3 than for unmethylated H3K423, indicating that H3K4me3-poised genes might be activated more efficiently or synchronously21. Detailed functional studies are needed to test these models and determine how the chromatin signature of embryonic pluripotency is established and regulates development.
ChIP was performed as described28 with minor modifications. For sequential ChIP, two changes were made to the protocol: (i) antibody was crosslinked to tosylactivated magnetic Dynabeads (Invitrogen) and (ii) chromatin was eluted from the beads in a milder elution buffer. After removal of the beads, 10μl of each sample was used for analysis of the first ChIP. The remaining 100μl was used for sequential ChIP. The eluate was diluted 10-fold in IP reconstitution buffer to reconstitute the original IP buffer composition. Samples were then pre-cleared by incubation with empty beads. After removal of the beads, half of each sample was subjected to a second round of immunoprecipitation with the tosylactivated beads crosslinked to the antibody against the ‘other’ modification (i.e. the one that was not used in the first ChIP) and the other half was incubated with empty beads as a negative control. Washing, elution and reversal of the crosslinks were performed as in single ChIP. See also Supplementary Figure 8. All sequential ChIPs were normalized to their input, i.e. the signal of the modification analyzed in the first round of ChIP, and fold-enrichment was calculated over a negative control.
For analysis on tiling arrays, ChIP material was amplified as previously described29 using the Sigma GenomePlex Kit to perform one round of amplification. Sample labeling, hybridization, washing and scanning were performed as described in the NimbleChip Arrays User’s Guide. Each array was normalized by MA2C13. Each probe was assigned a MA2Cscore to reflect the normalized and window-averaged log2 ratio of ChIP enrichment over control. See Supplementary Figure 1 for decision trees used to determine the chromatin and RNA pol II status of genes.
Zebrafish were maintained and raised under standard conditions. Wild-type and UAS14:GFP31 embryos were collected at the 1-cell stage, synchronized and allowed to develop to the desired stage at 28°C. GAL4-VP16 RNA was injected at the 1-cell stage at 48 pg per embryo.
For Western blotting, embryos were dechorionated immediately after fertilization, allowed to develop to the desired stage, deyolked either manually or, in case of H3K36me3 mechanically32, and snap frozen in liquid nitrogen. For every antibody, equal numbers of embryos were loaded for each developmental stage (RNA pol II (n=2), H3K36me3 (n=30), H3K4me3 (n=5) and H3 (n=1)). Samples were run on 6% (RNA pol II) or 15% (histones) acrylamide gels. Tubulin was analyzed on all blots as a loading control. Western blots with antibodies that recognize phosphoserine 5 or 2 on the C-terminal domain of RNA pol II indicated that the bandshift in Figure 1a is due to phosphorylation of the C-terminal domain (data not shown).
Antibodies used were Histone H3 (abcam ab1791), H3K4me3 (abcam ab8580 and millipore #07-473), H3K27me3 (upstate #07-449), H3K36me3 (abcam ab9050), RNA pol II (8WG16/Covance MMS-126R), phosphoserine 5 version of RNA pol II (H14/abcam ab24759), phosphoserine 2 version of RNA pol II (H5/abcam ab24758), and α-tubulin (Sigma T6074).
ChIP was performed essentially as described28 with some minor modifications. Briefly, for each condition, ~200 embryos were carefully staged, dechorionated and fixed in 1.85% formaldehyde for 15 minutes at room temperature.
Formaldehyde was quenched by adding glycine to a final concentration of 0.125M. Embryos were rinsed 3 times in ice-cold PBS, immediately resuspended in cell lysis buffer (10mM Tris-HCl pH7.5/10mM NaCl/0.5%NP40) and lysed for 15 min on ice. Nuclei were collected by centrifugation, resuspended in nuclei lysis buffer (50mM Tris-HCl pH7.5/10mM EDTA/1%SDS) and lysed for 10 min on ice. Samples were diluted 3 times in IP dilution buffer (16.7mM Tris-HCl pH7.5/167mM NaCl/1.2mM EDTA/0.01%SDS) and sonicated to obtain fragments of ~500 bp. Triton X-100 was added to a final concentration of 0.75% and the lysate was incubated overnight while rotating at 4°C with 25 to 50μl of protein G magnetic Dynabeads (Invitrogen) that had been pre-bound to an excess amount of antibody. Bound complexes were extensively washed with RIPA (50mM HEPES pH7.6/1mM EDTA/0.7%DOC/1%Igepal/0.5MLiCl) and TBS and then eluted from the beads with elution buffer (50mM NaHCO3/1%SDS). Cross-links were reversed o/n at 65°C, RNA degraded by a 30 min RNAseA treatment and DNA purified by QIAquick PCR purification kit (QIAGEN).
The first round of ChIP was essentially as described above, except for two changes: (i) Tosylactivated magnetic Dynabeads (Invitrogen) were used to which antibody was crosslinked according to the protocol provided by the manufacturer and (ii) chromatin was eluted from the beads in a milder elution buffer (50mM NaHCO3/0.1%SDS). Elution was performed at RT for 15 min in 120μl elution buffer. After removal of the beads, 10μl was removed from each sample for subsequent analysis of the first ChIP. The volume was adjusted to 300μl with standard elution buffer (1% SDS), and this sample was then used in the single ChIP protocol. The remaining 100μl was used for sequential ChIP. The eluate was diluted 10-fold in IP reconstitution buffer (30.9mM Tris pH7.5/123mM NaCl/4.6mM EDTA/0.37% SDS) and Triton X-100 was added to a final concentration of 0.75% to reconstitute the original IP buffer composition. Samples were then pre-cleared by incubating them with empty beads for 30 min at 4°C. After removal of the beads, half of each sample was subjected to a second round of immunoprecipitation with the tosylactivated beads crosslinked to the antibody against the ‘other’ modification (i.e. the one that was not used in the first ChIP) and the other half was incubated with empty beads as a negative control. Washing, elution and reversal of the crosslinks were performed as in single ChIP. See also Supplementary Figure 8. All sequential ChIPs were normalized to their input, i.e. the signal of the modification analyzed in the first round of ChIP. Fold-enrichment was calculated over a negative control; hoxa13b-ups served as a negative control for H3K27me3 > H3K4me3 sequential ChIPs, because this hoxa13b upstream region is occupied by H3K27me3 but not H3K4me3, and β-actin 1 served as a negative control for H3K4me3 > H3K27me3 sequential ChIPs, because this region is occupied by H3K4me3 but not H3K27me3. See Supplementary Discussion 2 for details.
For analysis on tiling arrays, ChIP material was amplified as previously described29 using the Sigma GenomePlex Kit to perform one round of amplification. Sample labeling, hybridization, washing and scanning were performed as described in the NimbleChip Arrays User’s Guide. H3K4me3 (abcam ab8580), H3K27me3 (upstate #07-449) and H3K36me3 (abcam ab9050) ChIPs were normalized to H3 (ab1791), RNA pol II (8WG16/Covance MMS-126R) ChIPs were normalized to input.
Sequences for each genomic region were retrieved from the Sanger Institute Zv7 sequence via the UCSC genome browser. Repetitive sequence elements were removed by RepeatMasker. 385K 50mer probes were printed using NimbleGen technology (NimbleGen, Madison, WI) and were spaced at 30 bp, resulting in 80 bp resolution. Tiled regions include specifically selected genes (i.e. maternally provided genes, genes that are transcribed at zygotic genome activation, genes important in development and genes known to be expressed in later stages) as well as two large contiguous genomic regions on Chromosomes 3 and 11. A comparison of the data obtained using all genes on the array with the data obtained using only the genes on the two contiguous regions resulted in very similar histone modification profiles showing that the array is a good representation of the genome. All coprecipitated DNA fragments from ChIP experiments were labeled and hybridized to the custom-designed array using the standard NimbleGen protocol (Madison, WI) as described12. Probe intensities were extracted using the standard ChIP analysis settings provided in the Nimblescan® software (NimbleGen, Madison WI).
After scanning, each array was normalized by MA2C13. Each probe was assigned a MA2Cscore to reflect the normalized and window-averaged log2 ratio of ChIP enrichment over control. One replicate was performed at 256-cell and the absence of signal was confirmed by ChIP-qPCR analysis at that stage (Supplementary Figure 7A). Two biological replicates were performed at dome/30% stage. The overall Pearson correlations for MA2Cscores of all probes between two biological replicates at dome/30% were 0.91 (H3K4me3), 0.92 (H3K27me3), 0.91 (H3K36me3), and 0.78 (RNA pol II), respectively. Since there may be minor staging differences between replicates, the data from one replicate (the first) was used in all subsequent analyses. Analysis of both replicates gave very similar results (Supplementary Figure 7B). Enriched histone modification/RNA pol II domains were determined with the following approach: i) MA2C was applied to detect enriched peaks under two different threshold cutoffs: FDR 5% and FDR 20%; i.e. each peak with MA2Cscore of x is assigned an FDR based on (the number of control peaks over ChIP with a score above x) divided by (the number of ChIP peaks over control with a score above x). ii) Adjacent peaks (under FDR 20%) were concatenated if the gap-size separating them was smaller than 5 kb and there were either no probes in the gap, or if at least 75% of the probes in the gap had a MA2Cscore larger than 0. iii) Domains that were smaller than 1.5 kb and did not contain peaks under FDR 5% were discarded.
A total of 822 RefSeq genes were used in the analysis (this included all the RefSeq genes larger than 1.5 kb that were tiled on the array from at least 1.5 kb upstream of the transcription start site up to the transcription termination site [−1.5kb to TSS, TTS]). Enriched chromatin domains (as defined above) were associated with these RefSeq genes. By following the decision trees shown in Supplementary Figure 1, this resulted in 3 groups of genes for RNA pol II, H3K36me3 and H3K4me3 (high, moderate, and low), and 4 groups of genes for H3K27me3 (high whole gene (W), high promoter (P), moderate, and low).
The metagene profiles were generated from the average MA2Cscores across all RefSeq genes (except in Supplementary Figure 2, where only genes with the indicated profile are included) by CEAS (Cis-regulatory Element Annotation System)33. In the profiles, transcription units are shown as metagenes, i.e. relative distance from TSS to TTS, whereas immediate upstream and downstream sequences are shown in absolute distance (bp).
Mean MA2Cscores were used to generate the heat map. For H3K36me3, the mean MA2Cscore of the last 2/3 of the gene body was used, while for H3K4me3, and H3K27me3, the mean MA2Cscores of the promoter regions ([−1.5kb, 1.5kb] to TSS) were used. To ensure unambiguous grouping of genes, genes with a moderate status in any histone modification (H3K4me3, H3K27me3, H3K36me3) were filtered out.
Scatter plots were generated to display the pairwise correlation between ChIP-chip data. For this, the mean MA2Cscore was used. For H3K36me3, the mean MA2Cscore of the last 2/3 of the gene body was used for each gene, while for the other data, the mean MA2Cscores of the promoter region ([−1.5k, 1.5k] to TSS) was used.
Expression analysis was performed using a NimbleGen expression array containing 392,778 probes with a maximum of 12 probes from 37,157 known and predicted zebrafish genes. Groups of 12 embryos from the indicated time-points were homogenized in Trizol and RNA was extracted using the standard Trizol protocol. 1 ug of the resulting total RNA was amplified into aRNA using the Ambion MessageAmpII aRNA kit. Fluorescent labeling of samples was performed by random-prime reverse transcription of 10 ug input aRNA with either Cy3- or Cy5- dNTPs. Arrays were hybridized and washed per standard NimbleGen protocol. Microarray scanning was performed on an Axon Genepix scanner. Expression microarray quality control and data analysis was performed using NimbleScan and R/Bioconductor with RMA (quantile) normalization. To determine if a gene was maternally provided, the averaged expression data from 1- and 4-cell stage embryos was used. If a gene’s average expression level at these stages was larger than the mean expression value of all genes, the gene was regarded as maternally provided, and removed from the expression analysis of dome/30% sample. A different cut-off that eliminated fewer genes (i.e. defined a maternal gene if its expression was larger than twice the mean expression value of all genes) was also tested. This analysis produced similar results, demonstrating that the more stringent cut-off we used in our analysis eliminates most or all maternally provided genes.
We thank members of the Schier lab for help and advice, H. G. Shin, L. Taing and Z. J. Wu for computational analysis and helpful discussions, N. Follmer and B. Lilley for technical advice, and J. Dubrulle, N. Francis, R. Losick, S. Mango, T. van Opijnen and W. Talbot for discussions and critical reading of the manuscript. This work was supported by NIH grants to XSL (1R01 HG004069) and AFS (5R01 GM56211), and by EMBO and HFSP (LT-00090/2007) fellowships to NLV.
Author information. ChIPchip data is available under GEO accession number GSE20023. Reprints and permissions information is available at www.nature.com/reprints. The authors declare no competing financial interests.
Author contributions. NLV and AFS designed the study. NLV performed the experiments. YZ performed computational analysis. NLV, YZ, JR, XSL and AFS designed and performed data analysis. IGW provided technical support. FI provided RNA profiling data. AR provided analytical advice. NLV and AFS interpreted the data and wrote the paper with support from co-authors.