extfieldsofstudy:
Materials Science
Physics

text:
Sequence design-based control of DNA droplets formed from phase separation of DNA nanostructures

DNA has the potential to realize a controllable liquid-liquid phase separation (LLPS) system, because the design of its base sequences results in programmable interactions. Here, we have developed a novel DNA-based LLPS system which enables us to create 'DNA droplets' and to control their dynamic behaviour by designing sequences of the DNA nanostructure. We were able to change the phase separation temperature required for the formation of DNA droplets by designing the sequences. In addition, the fusion, fission, and formation of Janus-shaped droplets were controlled by sequence design and enzymatic reactions. Furthermore, modifications of proteins with sequence-designed DNAs allowed for their capture into specific droplets. Overall, our results provide a new platform for designing the phase behaviour of macromolecular structures, and paves the way for new applications of sequence-designed DNA in the creation of cell-mimicries, synthetic membraneless organelles, and artificial molecular systems.

Liquid-liquid phase separation (LLPS) of water-soluble molecules induces the formation of macromolecular droplets in water. In polymer chemistry and soft matter physics, LLPS has been studied to elucidate underlying mechanisms of phase separation [1][2][3] and to construct cell-like structures 4 or bioreactors 5 . With recent advances in LLPS in biological research fields, studies have found that macromolecular droplets comprised of proteins and RNAs are formed within living cells via moderate interactions among the biomolecules [6][7][8][9][10][11] . LLPS of biomolecules has garnered considerable attention as the droplets have an important role in regulating several biological functions 12 . Although several in vitro studies have reproduced droplet formation by changing the concentrations of proteins/RNAs and salts 6,9,13,14 , no study has demonstrated the regulation of droplet formation depending on the biological information encoded in the biopolymer sequence. The sequence of biopolymers is an essential factor influencing droplet formation 15 because it determines intra-and inter-molecular interaction strength 16 . However, the potential of regulating LLPS by the sequence design of biopolymers remains unexplored, although it is expected that such a technology would provide new insights into not only LLPS mechanisms, but also LLPS applications.
Considerable progress is being made in modulating molecular interactions by biopolymer sequence design in the field of DNA nanotechnology, where sequence-specific hybridization of nucleic acids plays an important role. DNA nanotechnology is a promising avenue for the construction of nanostructures [17][18][19][20][21][22][23][24] and for the investigation of their phase behaviour 25-27 because the self-assembly of DNAs is precisely predictable from sequence information. The nanostructures can be further assembled into micro-sized structures. For example, a DNA hydrogel can be formed from the selfassembly of branched-DNA nanostructures [28][29][30][31][32] , and it can serve as a substrate for cell-free protein expression 29 , cytoskeleton for artificial cells 30 , or carrier for drug delivery 33 . On the other hand, a few studies have reported the phase behaviour of DNA nanostructures, such as the phase separation of DNA nanostructures into DNA-poor and rich phases 25 , formation of a coacervate-like DNA liquid 26 , and re-entrant phase dynamics 27 . These studies focused on the effects of temperature, DNA concentration, and ionic strength. However, sequence dependency in the interaction between DNAs 3 has not been translated into the ability to regulate phase behaviour. Sequence-design-based phase behaviour of DNA nanostructures will enable us to develop a novel LLPS system with controllability at the molecular level.
We herein report the formation and control of a 'DNA droplet' made of sequence-designed DNAs ( Fig. 1). We first show that DNA nanostructures (Fig. 1a) exhibit the LLPS that results in droplet formation depending on the surrounding temperature (Fig. 1b). Next, we show how the sequence design influences the phase behaviour of DNA nanostructures. Then, we demonstrate that control of the dynamic behaviour of DNA droplets, including fusion, fission, and segregation of DNA droplets can be achieved by rational design of DNA sequences (Figs. 1c-1e). Finally, as a potential application of our DNA-based LLPS system, we demonstrate the capture and partitioning of protein-cargos in DNA droplets in a sequence-specific manner (Fig. 1f). 4 5 Fig. 1. Droplet formation of DNA nanostructures based on phase separation. a, Schematic of the Y-shaped DNA nanostructure (Y-motif), and strand sequences of the Y-motif. b, Temperaturedependent phase-changes of the Y-motifs. The Y-motifs are dispersed in solution at a high temperature ('dispersed phase'). Decreasing the temperature induces phase separation into a 'droplet phase' that is a liquid-liquid biphasic system composed of a Y-motif-rich liquid droplets phase in a continuous Y-motif-poor phase. By further decreasing the solution temperature, the phase changes from the droplet to a 'gel phase' that is a gel-liquid biphasic system composed of Y-motif-rich gel phase in the Y-motif-poor phase. The phase-change temperature is dominated by the sequence of DNAs. c, d, e, and f, Schematic representations for the control of DNA droplet and selective cargo capture and partitioning.

Formation of DNA droplets via the phase separation of DNA nanostructures
We designed a Y-shaped DNA nanostructure (hereafter referred to as Y-motif) with three 'sticky ends' (SEs) (Fig. 1a). The Y-motif is composed of three single-stranded DNAs (ssDNAs): Y1, Y2, and Y3 (Fig. 1a). The 'stems' of the Y-motif have 16 base pairs that hybridized at 75 ºC in a buffer condition (20 mM Tris-HCl (pH: 8.0), 350 mM NaCl) ( Supplementary Fig. 1); the SE was a palindromic 8 nucleotide (nt) sequence and the same for three branches of the Y-motif, allowing for connection between the Y-motifs.
We first observed the Y-motif solution using a confocal laser scanning microscope and a stage heater. The Y-motifs were stained with a dye molecule (SYBR Gold). At 70 °C, the Y-motifs were dispersed in a buffer solution ('dispersed phase') ( Fig. 1b) and no micrometre-sized structures were observed owing to dehybridization of the SE at this temperature ( Fig. 2a(I)). When the temperature was lowered to 62 °C, spherical microparticles (from several to several tens of micro-meters) were observed in the solution, although the temperature was still much higher than the melting temperature (Tm) of the 8 nt SE (Tm: 45.9 °C) ( Fig. 2a(II)). Interestingly, fusion of the two particles occurred when they collided with each other (Figs. 2b and 2d, and Supplementary Movie 1); the fusion behaviour is similar to that of the droplets formed via LLPS in living cells 8 . We called the particles that exhibit such fusion 'DNA droplets', where Y-motifs spontaneously phase separate into Y-motif-rich liquid droplets in equilibrium with a Y-motif poor phase 34 ('droplet phase') (Fig. 1b). The DNA droplets grew in size via collision with one another (Supplementary Fig. 2).
To evaluate the mobility of the Y-motifs in the droplets, fluorescence recovery after photobleaching (FRAP) experiments were performed. In the FRAP experiments, the Y-motifs were labelled using FAM-modified Y2 without SE to avoid the detection of the reattaching the dye molecules to the Y-motifs. The half region of droplets was photo-bleached and the changes in fluorescence intensity were measured in the bleached regions. The fluorescence intensity in the droplets was immediately recovered after the bleaching (Figs. 2f and 2h), which suggests that the Y-7 motifs can diffuse inside droplets (do not form static networks) and the dynamic dissociation and association of the Y-motifs result in the formation of the droplet.
When the temperature was further decreased to 25 °C ( Fig. 2a(III)), the microparticles did not fuse (Figs. 2c and 2e). In addition, in the FRAP experiments, the fluorescence intensity did not immediately recover at around this temperature (Figs. 2g and 2i). Thus, the Y-motifs formed static networks and the DNA droplets became Y-motif-rich hydrogels ('gel phase') ( Fig. 1b). A previous study using a microrheological technique reported that DNA nanostructures exhibited a gradual transition into a hydrogel with decreasing temperatures 27 . Our results suggest that DNA nanostructures transition through a droplet phase before gelation.
We further found that this phase behaviour between the dispersed, droplet, and gel phases was reversible. The hydrogels changed into droplets when temperatures increased to 62 °C ( Fig. 2a (IV)), and the Y-motifs dispersed at 70 °C ( Fig. 2a(V)). These temperature-dependent phase-changes were confirmed with further continuous temperature changes ( Supplementary Fig. 3).

Effects of the sequence and the number of branches on phase separations and changes
We next investigated the effects of the SE sequences on the phase-change temperatures between the dispersed and droplet phases (Td), and between the droplet and gel phases (Tg). We defined Td as the temperature at which droplet formation was observed ( Supplementary Fig. 4), and Tg as the temperature at which less recovery of fluorescence in the FRAP experiment was detected ( Supplementary Fig. 5). Td and Tg were identified by decreasing the temperature from 85 ºC at a rate of -1 ºC/min. The phase diagram (Fig. 3a)  for the Y-motifs whose SEs that were 4-12 nt long ( Supplementary Fig. 6), but not with those that were 2 nt SE (ΔH: -9.6 kcal/mol) ( Supplementary Fig. 7).
We also investigated the effects of the number of branches in the motifs on the Td and Tg. We additionally designed two motifs that had four or six branches with 8 nts long SEs (Supplementary Table 2 The fact that the DNA droplets were formed at higher temperatures than the Tm of the SEs suggests that unstable hybridization between the SEs plays a key role in the formation of the DNA droplets. Such weak interactions of the SEs would allow for the dynamic association/dissociation of the motifs (Figs. 2f and 2h), which leads to the fluidic behaviour of DNA droplets, e.g. fusion of the droplets.
On the other hand, the increase of the number of branches in a motif resulted in the increase of Td We considered that Tg was attributed to the stability of SEs and entanglement of motifs. Previous studies have reported that the gelation process in the branched-DNA nanostructures occurs at around the Tm of SEs 32 . On the other hand, our results show that Tg for the Y-motif with 4 nt SE was higher than the Tm of the 4 nt SE (Fig. 3a). DNA nanostructures can be concentrated by the phase separation, also decrease the mobility of the motifs. The Tg in motifs having four or six branches was higher than their Tm (Fig. 3c), and the Tg in the six SE motifs was higher than that in four. Therefore, the effects of the entanglement on Tg would become more dominant than the stability of the hybridization of the SEs.

Orthogonality in DNA droplets
We expected that the rational design of the motif sequence allowed for the control of the dropletdroplet interactions. We designed an 'orthogonal Y-motif' with 8 nt SE ( orth Y-motif) whose sequences were not complementary to the Y-motif with 8 nt SE, but had almost similar thermodynamic parameters (ΔH and Tm) in their stem and SEs (Supplementary Tables 3 and 4). The orthogonality of the SE sequences between the Y-and orth Y-motifs enables us to acknowledge the sequence-specific fusion of DNA droplets (Fig. 4a). The Y-and orth Y-motifs were labelled with FAM (green) and Alexa 405 (blue), respectively. Then all strands for both of the motifs were mixed in a solution, and the motifs were annealed on the stage heater from 85 °C at a rate of -1 °C/min to 60 °C.
As we expected, the Y-and orth Y-motifs were individually phase separated into droplets. Fusion events were observed between the droplets composed of the same motifs but not between the droplets  Table 5). The S-motifs play a role for the cross-bridging between the Y-and orth Ymotifs and can eliminate the orthogonality between the two types of motifs, which results in the formation of DNA droplets composed of Y-, orth Y-and S-motifs (Fig. 4d). The Cy5-labelled S-motif, FAM-labelled Y-motif, and Alexa405-labelled orth Y-motif were mixed in a 1:3:3 molar ratio and visualized at 65 °C. In the experiment, the fluorescence of Cy5, FAM, and Alexa405 were observed in one droplet, suggesting that the Y-and orth Y-motifs were mixed via the S-motifs in the droplets (Fig. 4e). The elimination of the orthogonality was also confirmed in a two-dimensional histogram of FAM/Alexa405 fluorescence intensities ( Supplementary Fig. 20). In addition, the droplets composed of Y-, orth Y-and S-motifs exhibited fusion (Fig. 4f) similar to the droplets of a sole motif composition. Thus, selective fusion of the DNA droplets was achieved in two motifs that were orthogonal and the addition of the motif for cross-bridging eliminated the selectivity.

Fission and segregation of DNA droplets
To further demonstrate the controllability of the DNA droplets, we designed a mechanism for their fission (Figs. 5a and 5b). DNA-RNA chimera strands were introduced into the S-motif, termed 'chimerized-S-motif (CS-motif)' that also has SEs for both the Y-and orth Y-motifs ( Fig. 5a and Supplementary Table 5). Two RNA regions were located at around the centre of the CS-motif and ribonuclease A (RNase A) decomposes the RNA parts. The enzymatic reaction results in splitting CS-motif into two three-branched portions with affinity to the Y-motif ('Y-portion') and the orth Ymotif (' orth Y -portion') ( Fig. 5a, Supplementary Fig. 21). After splitting, the cross-bridge ability of the CS-motif disappears, and the Y-and orth Y-motifs cannot be mixed owing to their intrinsic orthogonality. As a result, the DNA droplets exhibit micro-phase separation in a spinodal decomposition manner, and finally fission into multiple droplets composed of Y-motif/Y-portion or orth Y-motif/ orth Y-portion (Fig. 5b).
The Y-, orth Y-, and CS-motifs formed droplets, representing the mixing of the Y-and orth Y-motifs via CS-motifs. (Fig. 5c) By using S-motifs (non-cleavable by RNase A) in addition to CS-motifs (cleavable by RNase A), we demonstrated the formation of complex-shaped droplets. When the ratio of the CS-motif was large but a few S-motifs existed, addition of RNase A did not induce fission but caused the segregation of the droplets owing to the immiscibility and the partial cross-bridging by a few S-motifs (Fig. 5f).
Experimental results showed that in the 90 % ratio of the CS-motifs, the DNA droplets segregated 16 and finally formed Janus-shaped droplets that have two portions for the Y-and orth Y-motifs in the droplets ( Fig. 5g and Supplementary Figs. 25 and 26). Fluorescence of Cy3 ( orth Y-portions) and Cy5 (Y-portions) were colocalized to Alexa405 ( orth Y-motif) and FAM (Y-motif), respectively. Furthermore, the 50 % ratio of the CS-motifs resulted in the patchy-like pattern formation in the DNA droplets ( Supplementary Figs. 27 and 28). These results showed the existence ratio of the S-and CSmotifs determined the shape of the DNA droplets after the addition of the RNase A, and indicated that the S-motifs played a role for not only the cross-bridging but also for a 'surfactant' between the Y-and orth Y-motifs.

Selective cargo capture and partitioning in DNA droplets
The combination of our DNA-based LLPS system and other molecules would provide an application basis and increase the functionality of the DNA droplets. We focused on the capture of cargo into the DNA droplets via DNA-modification. DNA-streptavidin conjugates that have SEs for Y-or orth Y-motifs were prepared (Fig. 6a) and the selective capture of streptavidin in the DNA droplets was visualized (Fig. 6b).
DNA-modified streptavidin, labelled with DyLight 549 fluorophores, was added during the visualization of DNA droplets composed of Y-or orth Y-motifs. The streptavidin was successfully localized in specific DNA droplets, depending on the modified DNA sequence (Figs. 6c and 6d). In addition, the accumulation processes of the streptavidin into DNA droplets were also visualized ( Supplementary Fig. 29).
Moreover, by adopting the technique to prepare the Janus-shaped DNA droplets (Fig. 5f), partitioning of the cargos in DNA droplets was achieved (Fig. 6e-6j). Streptavidin, modified with SE sequences, were homogeneously distributed in DNA droplets in which the Y-and orth Y-motifs were mixed by the S-and CS-motifs ( Fig. 6f and 6i). After the addition of RNase A, DNA droplets segregated into the Janus shape, and streptavidin were partitioned into the complementary sides depending on the modified sequence ( Fig. 6g and 6j). As we demonstrated here, modification of targets with SE sequences, which encode interaction information, allowed us to capture them in DNA droplets and to control their partitioning in a sequence-dependent manner. These results suggest that DNA droplets have a potential to be carriers of designated molecules and to provide reaction environments for biochemical processes as artificial membraneless organelles.

Discussion
In this paper, we reported the sequence-dependent phase behaviour of the DNA nanostructures and the control of DNA droplets according to the information encoded in the SE sequences. The sequencedependent phase-changes were reversible between the dispersed, droplet, and gel phases in response to the surrounding temperature (Fig. 2a), which is similar to the gas-liquid-solid phase transition of matter. We expect that these findings will be applicable to other nanostructures such as DNA origami 22 , DNA-modified proteins 37 , or nanoparticles covered with DNA strands [38][39][40] . In other words, the phase behaviour shown in the present study could be expanded to matters that can be modified with DNAs, which may enable us to design the phase and create droplets for any material.
A remarkable point in our DNA droplet system is the controllability of the interaction and dynamic behaviour of the droplets (Figs. [3][4][5]. As living systems are well organized-dynamic structures whose behaviour is regulated by the information encoded in biopolymers 41 , our DNA-based LLPS system would provide a new basis for the development of artificial cell engineering. Additionally, our demonstration for the capture of proteins in DNA droplets with sequence specificity (Fig. 6) can also make significant contributions to synthetic biology to generate a biochemical reaction environment composed of designed molecules.
From the viewpoint of material science, the successful creation of the Janus-shaped DNA droplets ( Fig. 5g) would give us a fundamental technique to construct water-soluble complex-shaped droplets in an aqueous phase, which is intrinsically different from the conventional Janus-droplets 42 based on the water/oil immiscibility. Furthermore, the Janus DNA droplets were created owing to a surfactantlike role of DNA nanostructures for the cross-bridge of the two types of Y-motif. It suggests that multi-affinity surfactants would be achieved by the modification of DNAs to designated molecules.
Finally, we envision that their combination with biochemical circuits that can be programmed to produce and decompose specific DNA strands [43][44][45] , will lead to the design and control of autonomous dynamic behaviour of DNA droplets as if living organisms. Ultimately, such autonomous behaviour of macromolecular structures could serve for the development of molecular robotic systems 46, 47 , 23 which become comparable with living cellular systems in the future. We believe that this study is an important step toward the development of a novel LLPS system whose dynamical behaviour can be controlled by information encoded in component molecules.

Oligonucleotide preparation
Sequences for the DNA and DNA-RNA chimera strands are provided in the Supplementary Tables 1-3, and 5. The oligonucleotides were purchased from Eurofins Genomics (Tokyo, Japan) and were of OPC grade purification, except for the fluorophore-labelled and chimerized strands, that were purified to an HPLC grade. The oligonucleotides were dissolved in ultrapure water (18 MΩ cm in resistance) at 100 µM concentrations and stored at -20 ºC until use.

Repetitive phase-change of Y-motifs
DNA strands for the Y-motifs with 8 nucleotide long sticky ends were mixed in a test tube with the buffer, 5 µM for each strand. SYBR Gold was also added at 1 × concentration to the test tube to stain the DNA. The test tube was heated at 85 ºC for 3 min; then gradually cooled to 25 °C at a rate of -1 °C/min, using a thermal cycler (Mastercycler nexus X2, Eppendorf, Hamburg, Germany). The  Supplementary Fig. 3.

Measurement of the phase-change temperatures and FRAP experiments
Each motif was labelled with dye-modified DNA without sticky ends, that was mixed at a 10 % molar ratio in the solution. (e.g., 5 µM for Y-1 and Y-3, 4.5 µM for Y-2, and 0.5 µM for Y2_0_FAM). The sample solution was poured in the chamber on the stage heater and incubated at 85 °C for 3 min.
Subsequently, the temperature was decreased to a designated temperature at a rate of -1 °C/min. To determine the temperature at which the DNA droplets were formed, the sample was visualized at the 1 °C interval. To evaluate the phase-change from the droplet phase to the gel phase, the fluidity of the motifs in the particles formed by the phase separation was analysed using FRAP experiments.  Supplementary Table 5. Sequences of six-junction motifs (S-motif) and the chimerized-S-motif

Measurement of Tm for the Y-motif
DNA strands for the Y-motifs without sticky ends were mixed in a test tube at 5 µM of each strand with a buffer, consisting of 20 mM Tris-HCl (pH 8.0) and 350 mM NaCl. SYBR Green I (Takara Bio, Kusatsu, Japan) was added at 1× concentration to the test tube to detect the formation of the double helix. The test tubes were set on a real-time PCR detection system (CFX Connect, Bio-Rad, Hercules, CA, USA). The sample was heated at 85 ºC for 3 min; then cooled to 5 °C at a rate of -1 °C/min.
Subsequently, the temperature was increased from 5 °C to 85 °C at the same rate. The fluorescence intensity was measured at 1 °C intervals during the decrease and increase of temperature ( Supplementary Fig. 1).

Concentration estimation
The concentrations of the Y-motifs in the droplets and gel phases were estimated by measuring the fluorescence intensity of the supernatant, after centrifugation, that precipitates the DNA droplets and hydrogels. First, the calibration curve for the intensity-based concentration values was prepared using Y-motifs without sticky ends (composed of Y-1_0, Y-2_0, Y-3_0 and 10 %(mol) Y-2_0_FAM). The Y-motif without a sticky end was annealed at 5 µM in a 20 µl volume in the PCR test tube with the buffer, using the thermal cycler from 85 °C to 4 °C, at a rate of -1 °C/min. The annealed Y-motifs without sticky ends were then diluted to 100, 10, and 1 nM for a 20 µl volume using the buffer, and the fluorescence intensity was measured using a spectrometer (FP-6500, JASCO, Hachioji, Japan).
Next, the Y-motifs with 4, 6, 8, 10, and 12 nucleotide long sticky ends were also annealed by the same procedure. After annealing, to precipitate the DNA droplets and hydrogels at the bottom of the test tubes, the PCR tube was placed in a larger test tube as a support during centrifugation ( Supplementary   Fig. 18 a). Then, the tubes were centrifuged at 2000 ´ g at 4 °C for 20 min. The precipitation of the phase separated Y-motifs were confirmed using a gel imager (Gel Doc EZ system, Bio-Rad, Hercules, CA, USA) ( Supplementary Fig. 18 b). The PCR tubes containing the precipitate were set on the thermal cycler. For measuring the concentrations in the droplet and gel phases, the temperatures were set at 1 °C lower than the Td for each length of the sticky end or 4 °C, respectively. The test tubes were incubated for 10 min at each temperature, and the supernatant was collected using a gel loading tip (BM Equipment, Tokyo, Japan). Ten micro-litres of supernatant was diluted to 200 µl using the buffer, and the fluorescence intensity was measured. The concentration of Y-motifs was estimated using the calibration curve.

Observations of orthogonality in the DNA droplets and its elimination
For the orthogonality experiments, DNA strands were mixed in a test tubes with the buffer at the concentrations shown in Table 6. This mixture was poured in the observation chamber and set on the stage heater. Then, the temperature was kept at 85 °C for 3 min, and decreased at a rate of -1 °C/min.
The fusion events of the DNA droplets were visualized at 60 °C using the confocal microscope at a scanning velocity of 7 s per frame.
To confirm the elimination of the orthogonality, the Y-, orth Y-, and S-motif were mixed in a test tube with the buffer at the concentrations shown in Table 7. The mixture was then visualized using the same method as described above. Supplementary

Fission of DNA droplets
To achieve fission of the DNA droplets, part of the DNA strands constituting the S-motif were replaced with DNA-RNA chimera strands, considering the motif geometry for the split into the Yand orth Y-portions. DNA and DNA-RNA chimera strands were mixed at the concentrations given in Supplementary

Segregation of DNA droplets
The Janus-shape and patchy-like patterns in the DNA droplets, were created by adjusting the concentrations of DNA-RNA chimera strands. For the Janus-shape and patchy-like pattern, 90 and 50 mol % of the DNA strands were replaced with the chimera strands, respectively (Supplementary Table 9 and 10). Sample solutions in the test tubes were annealed from 85 °C to 25 °C at a rate of -1 °C/min using the thermal cycler; then, 1 µl of the RNase A solution was added to the solutions for the 20 µg/ml final concentration. These mixtures were poured into the observation chamber on the stage heater and the temperature was increased to 62 °C, and visualized using the confocal microscope. Supplementary

Supplementary Discussion
The phase diagram (Fig. 3a in the main text) showed that the phase-change temperature between the dispersed and droplet phase (Td) for each sticky end (SE) was higher than the Tm of the SEs. The Tm is generally defined as the temperature at which half of the DNA strands in a solution formed double helixes when two complementary DNA strands were present. On the other hand, our motifs have multiple SEs (e.g., Y-motif has three SEs). Thus, the Tm for the SEs should not be considered, but instead, the Tm for the motifs (Tm-Motif) at which half of the motifs in a solution are connected, should be considered.
In the branched motifs, the hybridization of one of multiple SEs would be sufficient for the connection between the motifs. Based on this idea, we considered that the Tm-Motif can be described as the temperature at which at least one of the multiple SEs can hybridize at a 50 % ratio. It is calculated by the following equation:  Fig. 13); this trend is consistent with the experimental results (Fig. 3c). However, there were still gaps between the Td and Tm-Motif, in particular, for the shorter SEs, which suggests that another parameter can influence the phase-change temperature.
The phase diagram also showed that the difference between the Td and Tg (the phase-change temperature between the droplet and gel phase) becomes larger with decreasing SE length. As the hybridization stability of the SEs would play a key role in determining the phase-change temperature of both the Td and Tg, the width of the melting range (Tw) provides an indication for the phase-change temperature.
In addition to Tm, Tw is also an important parameter for the characterization of melting behaviour.
Tw can be defined as the difference between the temperatures at the end and the start of DNA melting ( Supplementary Fig. 14). DNA strands transfer from ssDNA to dsDNA states in the range of Tw and vice-versa. At the start of the melting, SEs can exhibit weak interactions, and at the end of melting, SEs form stable double helixes. As the self-assembly of motifs into droplets and gel structures results from weak interactions and stable network-formations, we supposed that temperature differences between Td and Tg (ΔTdroplet-gel) related to Tw.
We calculated Tw for each of the SEs from the melting curve using NUPACK 1 (Supplementary Fig.   15). The ΔTdroplet-gel and Tw were shown in Supplementary Fig. 16. ΔTdroplet-gel and Tw exhibited closer values and the differences between ΔTdroplet-gel and Tw are smaller in the shorter SEs. Although the entanglement of motifs can influence Tg, as we discussed in the main text, this graph suggests that Tw may be one of the parameters that determines the phase-change temperature.
Taken together, the entanglements of motifs, multiple SEs in motifs that provide new parameters of Tm-Motif, and Tw of SEs could be possible factors influencing the phase-change temperature and their combinatorial effects may determine Td and Tg. Model-based numerical simulations could be a


extfieldsofstudy: Biology Medicine

text:
Genome-Wide Identification and Expression Analyses of the bZIP Transcription Factor Genes in moso bamboo (Phyllostachys edulis)

The basic leucine zipper (bZIP) transcription factor (TF) family is one of the largest gene families, and play crucial roles in many processes, including stress responses, hormone effects. The TF family also participates in plant growth and development. However, limited information is available for these genes in moso bamboo (Phyllostachys edulis), one of the most important non-timber forest products in the world. In the present study, 154 putative PhebZIP genes were identified in the moso bamboo genome. The phylogenetic analyses indicate that the PhebZIP gene proteins classify into 9 subfamilies and the gene structures and conserved motifs that analyses identified among all PhebZIP proteins suggested a high group-specificity. Microsynteny and evolutionary patterns analyses of the non-synonymous (Ka) and synonymous (Ks) substitution rates and their ratios indicated that paralogous pairs of PhebZIP genes in moso bamboo underwent a large-scale genome duplication event that occurred 7–15 million years ago (MYA). According to promoter sequence analysis, we further selected 18 genes which contain the higher number of cis-regulatory elements for expression analysis. The result showed that these genes are extensively involved in GA-, ABA- and MeJA-responses, with possibly different mechanisms. The tissue-specific expression profiles of PhebZIP genes in five plant tissues/organs/developmental stages suggested that these genes are involved in moso bamboo organ development, especially seed development. Subcellular localization and transactivation activity analysis showed that PhebZIP47 and PhebZIP126 were localized in the nucleus and PhebZIP47 with no transcriptional activation in yeast. Our research provides a comprehensive understanding of PhebZIP genes and may aid in the selection of appropriate candidate genes for further cloning and functional analysis in moso bamboo growth and development, and improve their resistance to stress during their life.

Introduction
Transcription factor (TF) proteins play crucial roles in activating gene expression by binding to the promoter regions of target genes. The basic leucine zipper (bZIP) domain encoding genes, generally associated with stress resistance and control many crucial biological processes in plants, form one of the largest TF families in eukaryotes, and are named for a conserved bZIP domain containing a basic region and a leucine zipper sequence of 60-80 amino acid residues [1]. The basic region has a conserved N-x 7 -R/K-x 9 motif, which is responsible for DNA binding and nuclear localization, while the leucine zipper region is composed of several heptad repeats of leucine residues or other bulky hydrophobic amino acid, such as isoleucine, valine, phenylalanine, or methionine [2].
Moso bamboo (Phyllostachys edulis) belongs to the Poaceae family, represents the only major lineage of grasses that is native to forests, is one of the most important non-timber forest products and has great ecological, cultural, and economic value in the world [34]. It is used as an important source of timber, food, and paper material in daily life [34]. Extreme climate and unfavorable growth conditions such as excess salinity, soil depletion, drought, cold and heat stress, and wounding largely limit the growth of moso bamboo. The moso bamboo has evolved various adaptive mechanisms to deal with varied, adverse environmental factors, and survive those stresses. To our knowledge, however, very limited information is available regarding the bZIP TF gene family in moso bamboo, its multiple roles in stress and hormone effects not yet studied, until now. We want to preliminarily explore the function of bZIP TF gene in the stress and hormone response, and whether these genes play an important role in regulating growth and development as well as stress resistance in moso bamboo during their life.
We identified 154 putative bZIP TF genes (PhebZIPs) from the P. edulis genome, and investigated genomic structure and phylogenetic relationships with counterparts from the other Poaceae species O. sativa and B. distachyon. We systematically analyzed conserved motifs, microsynteny and evolutionary patterns of the encoded PhebZIP proteins, and located cis-regulatory-elements in the genomic sequences. According to the number of the cis-regulatory elements, 18 genes were selected to examine the expression of PhebZIP genes in response to three different phytohormones, as well as in different tissues/organs/developmental stages. For deeper study, we cloned two genes from 18 genes and analyzed their subcellular localization. In addition, the PhebZIP47 gene of these two genes was highly expressed in all hormone treatments and tissues was selected for transactivation activity analysis. In general, we provide a basis for the screening of appropriate candidate PhebZIP genes for further cloning to know the molecular mechanisms in hormone response and improve the stress resistance of moso bamboo.

Identification of bZIP Genes in moso bamboo
The candidate PhebZIPs were named based on chromosomal locations, PhebZIP1 to PhebZIP154, from top to bottom, in the present study. The proteins encoded by the 154 predicted full-length PhebZIP genes range from 80 (PhebZIP15) to 839 (PhebZIP134) amino acids (aa), and from 8679.74 to 91904.76 Da in relative molecular weight. The isoelectric points (PIs) of the PhebZIP proteins are predicted to range from 4.7 (PhebZIP28) to 11.73 (PhebZIP141). Detailed information for the examined bZIP proteins is shown in Table S1.

Phylogenetic Analysis of bZIP Proteins in the Three Different Plant Species
A comprehensive phylogenetic tree containing 339 bZIP gene protein sequences from three monocot species (96 from B. distachyon, 89 from the O. sativa, and 154 from P. edulis) was estimated ( Figure 1). The OsbZIP and BdbZIP sequences were obtained, and we generated a phylogenetic tree using the Neighbor-joining (NJ) algorithm in MEGA software (MEGA6.0) with the PhebZIP sequences. Detailed bZIP TF gene information from rice and B. distachyon is given in Table S2. According to the comprehensive phylogenetic tree, the result showed that the predicted PhebZIP gene family cluster into 9 subfamilies, named groups B-J. The three tested monocot species had representatives in most all subfamilies, group A only contains B. distachyon sequence. predicted to range from 4.7 (PhebZIP28) to 11.73 (PhebZIP141). Detailed information for the examined bZIP proteins is shown in Table S1.

Phylogenetic Analysis of bZIP Proteins in the Three Different Plant Species
A comprehensive phylogenetic tree containing 339 bZIP gene protein sequences from three monocot species (96 from B. distachyon, 89 from the O. sativa, and 154 from P. edulis) was estimated ( Figure 1). The OsbZIP and BdbZIP sequences were obtained, and we generated a phylogenetic tree using the Neighbor-joining (NJ) algorithm in MEGA software (MEGA6.0) with the PhebZIP sequences. Detailed bZIP TF gene information from rice and B. distachyon is given in Table S2. According to the comprehensive phylogenetic tree, the result showed that the predicted PhebZIP gene family cluster into 9 subfamilies, named groups B-J. The three tested monocot species had representatives in most all subfamilies, group A only contains B. distachyon sequence. Phylogenetic tree and distribution of bZIP protein from three plant species: rice, Brachypodium distachyon, and moso bamboo. The phylogenetic tree was constructed using the neighbor-joining method as implemented in MEGA6.0 from a bZIP protein sequence alignment. Bootstrap values from 1000 replicates are displayed at each node. The proteins on the tree can be divided into ten clades and the different clades are indicated by different colors. Phylogenetic tree and distribution of bZIP protein from three plant species: rice, Brachypodium distachyon, and moso bamboo. The phylogenetic tree was constructed using the neighbor-joining method as implemented in MEGA6.0 from a bZIP protein sequence alignment. Bootstrap values from 1000 replicates are displayed at each node. The proteins on the tree can be divided into ten clades and the different clades are indicated by different colors.

PhebZIP Gene Structure, Conserved Protein Motifs
The exon-intron structure of the 154 PhebZIP genes was examined to obtain insight into the family's structural evolution ( Figure 2). The results showed that the number of introns of PhebZIP genes varied from 0 to 17. Interestingly, no intron was detected in PhebZIP genes of subfamily I and H. In summary, most of PhebZIPs belonging to the same groups showed similar exon-intron structures, which indicates their close evolutionary relationship and the classification of subfamilies. The exon-intron structure of the 154 PhebZIP genes was examined to obtain insight into the family's structural evolution ( Figure 2). The results showed that the number of introns of PhebZIP genes varied from 0 to 17. Interestingly, no intron was detected in PhebZIP genes of subfamily I and H. In summary, most of PhebZIPs belonging to the same groups showed similar exon-intron structures, which indicates their close evolutionary relationship and the classification of subfamilies.
The MEME Web server was used to investigate the diversity of motifs in the translations from our candidate PhebZIP genes. MEME motifs were captured and designated motif 1 to motif 20 ( Figure 3). Notably, almost all of the PhebZIPs translations contain motif 1 and subfamily F only contains motif 1. Furthermore, motif 17 is only found in subfamily B and motif 8 also presents exclusively in subfamily G, indicating that these motifs might have a special function. In general, many of the PhebZIPs encoded proteins of the same subfamily contain very similar motif compositions and they might have had similar functions.  The MEME Web server was used to investigate the diversity of motifs in the translations from our candidate PhebZIP genes. MEME motifs were captured and designated motif 1 to motif 20 ( Figure 3). Notably, almost all of the PhebZIPs translations contain motif 1 and subfamily F only contains motif 1. Furthermore, motif 17 is only found in subfamily B and motif 8 also presents exclusively in subfamily G, indicating that these motifs might have a special function. In general, many of the PhebZIPs encoded proteins of the same subfamily contain very similar motif compositions and they might have had similar functions.

Chromosomal Location, Microsynteny and Evolutionary Patterns of PhebZIP Genes
Gene duplication events play a vital role in the evolution of genes and can lead to the appearance of orthologous and paralogous pairs. According to the chromosomal location map ( Figure S1), 154 bZIP genes of moso bamboo were distributed irregularly across 23 chromosomes. To further investigate the evolution of the PhebZIP gene family, using microsynteny analysis, we investigated the genome duplication events within moso bamboo ( Figure 4). A total of 132 collinear pairs were identified in moso bamboo (Table S3). To estimate evolutionary rates and determine the relative divergence of moso bamboo, the Ks values and Ka/Ks ratios were estimated for the paralogous (Pe-Pe) gene pairs ( Figure 5). The frequency distribution of relative Ks values peaked at 0.1-0.2 ( Figure 5A), suggesting that a large-scale genome duplication event occurred around 7-15 MYA creating the moso bamboo bZIP TF gene paralogous pairs. A Ka/Ks value of < 1 indicates that a gene underwent negative or purifying selection, while Ka/Ks =1 and > 1 indicate neutral selection, and positive selection, respectively. All the bZIP TF gene paralogue Ka/Ks ratios in moso bamboo are < 1 ( Figure 5B), indicating that the PhebZIP family paralogous principally underwent purifying selection.

Chromosomal Location, Microsynteny and Evolutionary Patterns of PhebZIP Genes
Gene duplication events play a vital role in the evolution of genes and can lead to the appearance of orthologous and paralogous pairs. According to the chromosomal location map ( Figure S1), 154 bZIP genes of moso bamboo were distributed irregularly across 23 chromosomes. To further investigate the evolution of the PhebZIP gene family, using microsynteny analysis, we investigated the genome duplication events within moso bamboo ( Figure 4). A total of 132 collinear pairs were identified in moso bamboo (Table S3). To estimate evolutionary rates and determine the relative divergence of moso bamboo, the Ks values and Ka/Ks ratios were estimated for the paralogous (Pe-Pe) gene pairs ( Figure 5). The frequency distribution of relative Ks values peaked at 0.1-0.2 ( Figure 5A), suggesting that a large-scale genome duplication event occurred around 7-15 MYA creating the moso bamboo bZIP TF gene paralogous pairs. A Ka/Ks value of < 1 indicates that a gene underwent negative or purifying selection, while Ka/Ks =1 and > 1 indicate neutral selection, and positive selection, respectively. All the bZIP TF gene paralogue Ka/Ks ratios in moso bamboo are < 1 ( Figure 5B), indicating that the PhebZIP family paralogous principally underwent purifying selection.

Analysis of PhebZIP cis-Regulatory Elements
Gene function and regulation are largely determined by cis-regulatory elements [35]. Three types of cis-acting elements involved in gibberellin, abscisic acid response and MeJA were identified in PhebZIP gene regulatory regions in our study (Table S4). The cis-regulatory element analysis shows that the 18 selected PhebZIP genes contain many elements associated with GA, ABA and MeJA treatments ( Figure 6). Regarding hormone-related cis-acting elements, gibberellin responsive elements (GARE-motif, P-box, and TATC-box), abscisic acid responsive elements (ABRE), and MeJA stress responsive elements (CGTCA-motif) are found in the promoters of 10, 17, and 15 PhebZIP genes, respectively [36][37][38]. Notably, PhebZIP88 possesses 13 ABRE elements, PhebZIP68 and PhebZIP142 have 11 ABRE elements, respectively. There are 14, 96 and 50 elements that were associated with GA, ABA, and MeJA, respectively. Since these genes had the higher number of cis-regulatory elements related to GA-, ABA-and MeJA-responses; were further selected for qRT-PCR analysis.

Analysis of PhebZIP cis-Regulatory Elements
Gene function and regulation are largely determined by cis-regulatory elements [35]. Three types of cis-acting elements involved in gibberellin, abscisic acid response and MeJA were identified in PhebZIP gene regulatory regions in our study (Table S4). The cis-regulatory element analysis shows that the 18 selected PhebZIP genes contain many elements associated with GA, ABA and MeJA treatments ( Figure 6). Regarding hormone-related cis-acting elements, gibberellin responsive elements (GARE-motif, P-box, and TATC-box), abscisic acid responsive elements (ABRE), and MeJA stress responsive elements (CGTCA-motif) are found in the promoters of 10, 17, and 15 PhebZIP genes, respectively [36][37][38]. Notably, PhebZIP88 possesses 13 ABRE elements, PhebZIP68 and PhebZIP142 have 11 ABRE elements, respectively. There are 14, 96 and 50 elements that were associated with GA, ABA, and MeJA, respectively. Since these genes had the higher number of cis-regulatory elements related to GA-, ABA-and MeJA-responses; were further selected for qRT-PCR analysis.   Based on functional annotations, the cis-acting elements classify into five major groups: gibberellin, abscisic acid, and MeJA.

PhebZIP Gene qRT-PCR Analysis
The cis-regulatory elements analysis indicates that the 18 selected PhebZIP genes contain a large number of elements which are related to ABA-, GA-and MeJA-responses. Therefore, these results led us to wonder whether these selected genes were up-regulated in these hormone treatments and in different tissues/organs/developmental stages by qRT-PCR analysis.
The MeJA application can stimulate the expression of defense plant genes. Under the treatment, the expression of PhebZIP21, −26, −47, −68, −97, −105 and −144 was up-regulated, while PhebZIP59 was down-regulated, compared with 0 h (Figure 7). Furthermore, the expression of PhebZIP47 and PhebZIP68 peaked at 1 h; PhebZIP21, and PhebZIP26 peaked at 3 h; and PhebZIP105 and PhebZIP144 peaked at 6 h. The peak expression of PhebZIP97 was observed at 12 h.

PhebZIP Gene qRT-PCR Analysis
The cis-regulatory elements analysis indicates that the 18 selected PhebZIP genes contain a large number of elements which are related to ABA-, GA-and MeJA-responses. Therefore, these results led us to wonder whether these selected genes were up-regulated in these hormone treatments and in different tissues/organs/developmental stages by qRT-PCR analysis.
The MeJA application can stimulate the expression of defense plant genes. Under the treatment, the expression of PhebZIP21, −26, −47, −68, −97, −105 and −144 was up-regulated, while PhebZIP59 was down-regulated, compared with 0 h (Figure 7). Furthermore, the expression of PhebZIP47 and PhebZIP68 peaked at 1 h; PhebZIP21, and PhebZIP26 peaked at 3 h; and PhebZIP105 and PhebZIP144 peaked at 6 h. The peak expression of PhebZIP97 was observed at 12 h.     (Figure 9). Additionally, PhebZIP36 showed the greatest up-regulation (about 300-fold). The paralogous gene pair also showed similar expression patterns in response to GA. Finally, one of these gene (PhebZIP47) is up-regulated under all three treatments.   Samples were obtained from young leaves, mature leaves, roots, stems, and seeds. Color scale erected vertically at the right side of the picture represents log10 expression values, green indicates lower and red higher transcript abundance. L1: mature leaves; L2: mature leaves; R: roots; St: stems; Se: seeds.

Subcellular Localization and Transactivation Activity Analysis
For deeper study, we cloned two genes (one was highly expressed in all hormone treatments and the other was only highly expressed under GA treatment) from 18 genes and analyzed their subcellular localization. In addition, the PhebZIP47 gene of these two genes was highly expressed in all hormone treatments and tissues was selected for transactivation activity analysis. As shown in Figure 11, the control treatment shows green fluorescence throughout the membrane and within the cell nucleus, and the PhebZIP47-GFP and PhebZIP126-GFP fusion genes' green fluorescence signals are specifically distributed in the nucleus, consistent with previous predictions that TF binding of specific cis-elements occurs in the nucleus.
All the plasmids were independently transformed into the Y2HGold yeast strain and exhibited visible white colonies on the SD/-Trp medium to explore the transcriptional activity of PhebZIP47. The results showed only the positive controls grew well and turned blue in the SD/-Ade/-His/-Trp/X-a-GAL medium. In contrast, the pGBKT7-PhebZIP47 and negative control group did not grow on this medium (Figure 12), indicating that this gene has no transcriptional activity. Samples were obtained from young leaves, mature leaves, roots, stems, and seeds. Color scale erected vertically at the right side of the picture represents log10 expression values, green indicates lower and red higher transcript abundance. L1: mature leaves; L2: mature leaves; R: roots; St: stems; Se: seeds.

Subcellular Localization and Transactivation Activity Analysis
For deeper study, we cloned two genes (one was highly expressed in all hormone treatments and the other was only highly expressed under GA treatment) from 18 genes and analyzed their subcellular localization. In addition, the PhebZIP47 gene of these two genes was highly expressed in all hormone treatments and tissues was selected for transactivation activity analysis. As shown in Figure 11, the control treatment shows green fluorescence throughout the membrane and within the cell nucleus, and the PhebZIP47-GFP and PhebZIP126-GFP fusion genes' green fluorescence signals are specifically distributed in the nucleus, consistent with previous predictions that TF binding of specific cis-elements occurs in the nucleus.
All the plasmids were independently transformed into the Y2HGold yeast strain and exhibited visible white colonies on the SD/-Trp medium to explore the transcriptional activity of PhebZIP47. The results showed only the positive controls grew well and turned blue in the SD/-Ade/-His/-Trp/X-a-GAL medium. In contrast, the pGBKT7-PhebZIP47 and negative control group did not grow on this medium (Figure 12), indicating that this gene has no transcriptional activity.

Discussion
Increasing evidence shows that the bZIP TFs play crucial roles in various developmental and physiological processes and involved in various plant hormone and stress responses in plants.
Much of the research into the functions of the bZIP genes have been characterized in model plants, especially in Arabidopsis thaliana [3], and some grass family (Gramineae) species including Oryza

Discussion
Increasing evidence shows that the bZIP TFs play crucial roles in various developmental and physiological processes and involved in various plant hormone and stress responses in plants.
Much of the research into the functions of the bZIP genes have been characterized in model plants, especially in Arabidopsis thaliana [3], and some grass family (Gramineae) species including Oryza

Discussion
Increasing evidence shows that the bZIP TFs play crucial roles in various developmental and physiological processes and involved in various plant hormone and stress responses in plants. Much of the research into the functions of the bZIP genes have been characterized in model plants, especially in Arabidopsis thaliana [3], and some grass family (Gramineae) species including Oryza sativa, Sorghum bicolor, Zea mays, Brachypodium distachyon, Triticum aestivum [4][5][6][7][8]. To our knowledge, however, very little is known about this family in moso bamboo, the only major lineage of grasses. Moso bamboo is one of the most important species for afforestation and timber and has great ecological, cultural, and economic value in the world [34]. Compared with other gramineous plants, research on the function and growth of bamboo has evolved slowly, especially in regard to its response to adverse environmental conditions and plant hormones. Therefore, we want to select a gene family and preliminarily screen some functional genes related to plant growth and development, and stress and hormone response to solve some problems that bamboo will encounter in its growth. The bZIP TF gene family, one of the largest TF families in the eukaryotes, is involved in many biological processes, including plant growth and development, and stress and hormone response. Therefore, we selected this gene family for further analysis.
In this report, we used phylogenetic analysis to classify bZIP TF genes into 10 groups-154 PhebZIP genes, along with 89 bZIP genes from rice, and 96 from Brachypodium distachyon (Figure 1). We also examined their evolutionary history and selected 18 genes have the higher cis-elements under GA, ABA and MeJA stress conditions. For deeper study, two genes were selected for subcellular localization and PhebZIP47 gene with high expression under three treatments and in five tissues were selected for transactivation activity analysis. The information above provides a comprehensive overview of the bZIP gene family in moso bamboo and lays the foundation for future functional characterization of individual PhebZIP genes.

The Distribution and Divergence Times of bZIP Gene Family
Previous research has found that intron numbers of rice and B. distachyon are 0-12 and 0-13, respectively [4,7]. And in our research, the intron/exon structural analysis showed that the number of introns in the PhebZIP genomic sequences varies from 0 to 17 ( Figure 2). This suggested that there were similar gene structure diversity of bZIP genes in different species. Although it has been reported that duplicated genes rarely diverge with respect to their biochemical function; gene duplication is considered to provide the raw material for evolution; duplicated genes may undergo substantial changes in their structures and/or regulatory mechanisms allowing them to assume novel roles [39,40]. The plant specificity is often the result selective gene loss or gain during evolution and can involve chromosomal rearrangements and fusions, such as gene duplication with the three principal types of exon/intron diversification: exon-intron gain/loss, exonization/pseudo-exonization and insertion/deletion [39]. Most of the PhebZIP genes in subfamilies I, G, E, F, D and B contain 0-6 introns; the rest subfamilies C, J and H, most of the PhebZIP genes have the most introns. The MEME conserved motif analysis shows that motif 17 and 18 are unique to subfamily D, motif 8 is only found in subfamily G, and subfamily F only contains motif 1. In addition, the motifs 1 is present in most of the PhebZIPs. Many of the PhebZIPs encoded proteins of the same subfamily contain very similar motif compositions.
According to the chromosomal location map ( Figure S1), 154 bZIP genes of moso bamboo were distributed irregularly across 23 chromosomes. To estimate evolutionary rates and determine the relative divergence of moso bamboo, the Ks values and Ka/Ks ratios were estimated for the paralogous (Pe-Pe) gene pairs ( Figure 5). The frequency distribution of relative Ks values peaked at 0.1-0.2 ( Figure 5A), suggesting that a large-scale genome duplication event occurred around 7-15 MYA creating the moso bamboo bZIP TF gene paralogous pairs. A Ka/Ks value of < 1 indicates that a gene underwent negative or purifying selection, while Ka/Ks =1 and > 1 indicate neutral selection, and positive selection, respectively [41,42]. All the bZIP TF gene paralogue Ka/Ks ratios in moso bamboo are < 1 (Figure 5B), indicating that the PhebZIP family paralogous principally underwent purifying selection.

PhebZIP Genes Contain Higher Number of cis-Elements Were Selected for qRT-PCR Analysis
Accumulated evidence has confirmed that many bZIP TFs participate in varied developmental processes, and stress and hormone responses, including abscisic acid (ABA) [22], gibberellin (GA) [23], ethylene [24], pathogen infection [25], drought [26], cold [27]. It is worth noting that gene function and regulation are largely determined by cis-regulatory elements [35]. The cis-elements analysis showed that the PhebZIP gene containing many promoter elements related to hormone and stress response. However, we ultimately selected 18 PhebZIP genes containing a large number of promoter elements associated with GA, ABA, and MeJA-response for further analysis by qRT-PCR.
From the cis-elements analysis, the gibberellin-responsive cis-regulatory elements (GARE motif, P-box, and TATC-box) are seen in 10 of the PhebZIP genes. Among these genes, nine genes have one or more elements, and are highly expressed under GA treatment. Only one gene (PhebZIP97) contains one element and are not highly expressed under GA treatment. On the contrary, one moso bamboo genes PhebZIP97 has no CGTCA-motif elements, but is highly expressed after MeJA treatment. This phenomenon suggests that the number of promoters is related to gene expression in respond to adverse conditions but not the only criterion, and the expression induction of genes are complex biological processes.
Moso bamboo will encounter all kinds of adverse environmental conditions during their life. MeJA as a plant hormone and signal molecule related to injury, is widely existed in plants [25]. Exogenous application can stimulate the expression of defense plant genes and induce chemical defense. Under MeJA, which simulated mechanical damage and insect feeding have similar effects, seven of PhebZIP genes were significantly up-regulated, suggested that these gene increase resistance to stress. The ABA plays a vital role in plants responding to abiotic stress. In Arabidopsis, ABF2, ABF3, and ABF4 were significantly up-regulated by ABA. Further functional analyses showed that these genes were involved in the regulation of ABA-mediated processes, such as root elongation and seed germination [43]. In our research, the relative expression showed that seven of selected PhebZIP genes to be significantly up-regulated, and one gene was down-regulated after ABA treatment, indicating the comprehensive response of bZIP genes to ABA in moso bamboo. GA, which are tetracyclic diterpenoid growth factors, are essential regulators in many aspects of plant development, including seed germination, stem elongation, and flowering [23]. The expression pattern analysis showed that almost all of the genes we selected were highly expressed under GA treatment. The results suggest that these genes may be largely involved in plant development and growth.
To further understand the expression of these genes in different tissues/organs, growth and development, five organs were selected for qRT-PCR and the result indicated that most of the genes are expressed in these tissues, possibly suggesting that these genes play important roles and have widespread involvement in plant growth and development, and is consistent with previous research in Arabidopsis [3]. Simultaneously, the expression pattern analysis of tissues also shows all PhebZIP genes to be highly expressed in seeds, which may play a vital role in moso bamboo seed maturation [18]. It is worth noting that the expression patterns of paralogous PhebZIP gene pairs under the three hormone treatments and in the five organs are very similar. This further proves that these paralogous gene pairs are functionally similar. In summary, the expression profiles of 18 selected PhebZIP genes under GA, ABA, and MeJA and in five organs revealed that these genes are widely involved in growth and development, stress and hormone response.

Subcellular Localization and Transactivation Activity Analysis
In this study, the subcellular localization of PhebZIP47 and PhebZIP126 in the nucleus is consistent with previous predictions of bZIP TF binding to cis-elements in the nucleus, thus, these two genes are typical nuclear localization transcription factors. It is worth noting that the PhebZIP47 was highly expressed under three hormone treatments and in all tissues. We want to know if this gene is a transcriptional activator. The transactivation activity analysis showed that PhebZIP47 had no transcriptional activation activity in yeast. Similar results were observed in rice that OsbZIP71 has no transcriptional activity in yeast [44]. This phenomenon may be explained by the activation of this protein possibly relying on some posttranslational modifications or needing to be activated by some unknown upstream proteins [45].

PhebZIP TF Database Searches
All of the annotated protein sequences from the moso bamboo genome database (http://server.ncgr. ac.cn/bamboo/index.php) were downloaded to identify genomic bZIP TFs [46]. Moreover, to avoid missing any moso bamboo bZIP genes, we used the Hidden Markov Model (HMM) of the bZIP domain (Pfam PF00170) to confirm them in our protein dataset using their Web-based HMM search tool [47,48]. The common domains of the Rice and B. distachyon gene families were used in Blastp (e-value = le −5 ) to identify corresponding moso bamboo gene families. Finally, all non-redundant candidate PhebZIP sequence translations were checked against the PFAM (https://pfam.xfam.org) (access in January 2019). databases to confirm the existence of the bZIP domain [49]. Redundant sequences and those without bZIP domains were manually removed from our dataset.

Microsynteny and Evolutionary Patterns of the bZIP Genes in moso bamboo
The duplication of bZIP genes were identified using MCScanX (http://chibba.pgml.uga.edu/ mcscan2/) (access in January 2019). [54]. First, whole-genome protein sequences from moso bamboo were searched using BLASTP with an E-value cutoff of 1e −10 , and identity > 50%. Then, the default parameters of MCScanX were used for detecting the synteny regions. Finally, based on information about collinear pairs and genetic location, circos was used to create collinear analysis diagrams [55]. The selective pressure and the occurrence of duplication events were estimated on the datasets by calculating non-synonymous (Ka) substitution and synonymous (Ks) rates between the collinear pairs. A synonymous substitution (Ks) is defined as a mutation in which a nucleotide base is replaced by a different base in a protein-coding region of a gene that does not result in an amino acid change in the encoded protein, while a non-synonymous substitution (Ka) results in a change in the amino acid sequence of a protein [56]. The non-synonymous and synonymous substitution rates were then calculated using DnaSP 5 to analyze gene duplication events [57,58]. The Ks value was used as a proxy to estimate the timing of large-scale duplication events [59], with Ks converted into a time framework using the formula T = Ks/2λ × 10 −6 Mya [60]. The clock-like synonymous substitution rate variable (λ) was specified at 6.5 × 10 −9 years, as per previous estimates [46].

Analysis of PhebZIP TF Putative Promoter Regions
We examined 2000-bp regions upstream of the transcriptional start in all of our candidate PhebZIP sequences to identify putative cis-elements in these promoter regions. We used the Web-based promoter search Plant-CARE (http://bioinformatics.psb.ugent.be/webtools/plantcare/html/) (access in January 2019). [61] to identify those putative cis-acting elements related to stress responses, hormone effects.

Plant Materials and Phytohormone Treatments
Moso bamboo seeds were identified and provided by the Guilin Forestry Bureau, Guilin, Guang Xi Province, China, and were grown in plastic containers in a greenhouse under a 14-h-light/10-h-dark photoperiod at 25 ± 2 • C and 60-80% humidity. For phytohormone treatments, the young leaves of three-month-old moso bamboo seedlings were sprayed with 10 µM ABA, 100 µM GA, and 100 µM MeJA, separately. All samples were collected at 0, 1, 3, 6, 12, and 24 h separately from seedlings after stress application and untreated leaves (0 h) were used as a control. To obtain tissue-specific expression profiles of PhebZIP genes, the tissues of young leaves, mature leaves, roots, stems, and seeds from seedlings were also collected. All plant samples were immediately frozen in liquid N 2 and stored at −80 • C for RNA extraction. Three biological and three technical replicates were performed to reduce experimental error.

RNA Extraction and qRT-PCR Analysis
Total moso bamboo leaf RNA under different stresses, and from different organs, was extracted with TRIzol reagent (Invitrogen, NH, USA) according to the manufacturer's instructions. The primer pairs for the 18 selected PhebZIP genes were designed using Primer Express 3.0, and the specificity of primers was checked using BLAST searches against a coding sequence (CDS) database downloaded from the moso bamboo genome database (http://server.ncgr.ac.cn/bamboo/index.php) (access in January 2019). [46]. The tonoplast intrinsic protein 41 (TIP41) gene, which is stably expressed in moso bamboo, was used as a reference to normalize expression data [62]. The qRT-PCR amplifications using SYBR Green Master Mix reagent (Applied Biosystems) were performed on an ABI 7300 Real-Time system (Applied Biosystems) in 20 µL reactions containing 1 µL of each gene-specific primer, 1 µL of cDNA sample, 7 µL ddH 2 O, and 10 µL SYBR Green Master Mix reagent (Applied Biosystems) [63]. All primers for amplification of PhebZIP genes are given in Table S5. The qRT-PCR amplification conditions were: 95 • C for 30 s, followed by 40 cycles of 95 • C for 10 s, 55 • C for 15 s, and 72 • C for 10 s. A melting curve analysis was performed to determine those primers with the highest specificity and efficiency. The relative expression levels were assessed based on the 2 −∆∆CT method [64].

Subcellular Localization Analysis of PhebZIP47 and PhebZIP126
The full-length PhebZIP47 and PhebZIP126 coding sequences were cloned and inserted into the pCAMBIAI1305 vector (TaKaRa, Beijing, China), which contains a GFP gene. The constructed PhebZIP47-GFP and PhebZIP126-GFP vectors were then separately transformed into Agrobacterium tumefaciens EHA105 (BIO-RAD, Hercules, CA, USA), and syringes were used to inject the suspensions into N. tabacum (tobacco) leaves. The PhebZIP47-GFP and PhebZIP126-GFP expressions were observed using confocal laser scanning microscopy (CarlZeiss LSM710).

Transactivation Activity Analysis of PhebZIP47
The transactivation activity analysis of PhebZIP47 protein in yeast was determined as per the previously described method [65]. DNA fragments containing the whole open reading frames of PhebZIP47 was cloned and then inserted into the pGBKT7 vector (Clontech, Palo Alto, CA, USA) to create the pGBKT7-PhebZIP47 construct. Then, the recombinant vector, the negative control pGBKT7 plasmids and the positive control pGBKT7-53+pGADT7-T were used to transform the AH109 yeast strain according to the manufacturer's protocol (Stratagene, Santa Clara, CA, USA). The transformed strains were grown on SD/-Trp and SD/-Trp/-His/-Ade/X-a-Gal plates at 30 • C for 3-5 days.

Ethics Approval and Consent to Participate
The seeds of moso bamboo were collected from Guilin in Guang Xi Province, China in June, 2018. And the seeds were provided and identified by the Guilin Forestry Bureau. All the materials of moso bamboo used and analyzed were public and available for non-commercial purposes. This article did not contain any studies with human participants or animals performed by any of the authors.

Conclusions
Our results represent the first genome-wide identification and analysis of the bZIP TF family in moso bamboo. This includes phylogenetic tree, protein motif prediction, gene structure, microsynteny and evolutionary patterns, all of which suggested a complex evolutionary history of this family in moso bamboo. The expression profile analyses of 18 selected PhebZIP genes demonstrate that these genes play pivotal roles in many processes, including stress response, hormone signaling, and in plant growth and development. This research provides a comprehensive understanding of PhebZIP genes and may aid in the selection of appropriate candidate genes for further cloning and functional analysis in moso bamboo growth and development and improve their resistance to biological and abiotic stress. Acknowledgments: Special thanks go for the anonymous reviewers for their enlightening comments in reformulating this paper.

Conflicts of Interest:
The authors declare no conflict of interest.


extfieldsofstudy: Medicine Chemistry

text:
Lovastatin increases the proliferation and osteoblastic differentiation of human gingiva-derived stem cells in three-dimensional cultures

Lovastatin is a cholesterol-lowering agent that also has effects of cell proliferation and apoptosis. The present study was performed to evaluate the effects of lovastatin on the proliferation and osteogenic differentiation of three-dimensional cell spheroids formed from human gingiva-derived stem cells (GDSCs) using concave microwells. GDSCs were plated on polydimethylsiloxane-based concave micromolds and grown in the presence of lovastatin at concentrations of 0, 2 and 6 µM. The morphology of the cells was viewed under an inverted microscope, and cell viability was determined with Cell Counting kit-8 on days 2, 7 and 14. Alkaline phosphatase activity assays were performed to evaluate the osteogenic differentiation on days 2 and 8. Alizarin red-S staining was also used to assess the mineralization of the stem cell spheroids at day 14. The results confirmed that GDSCs formed spheroids in concave microwells. No significant changes were noted with longer incubation time, and no significant differences in cell viability were noted between the three lovastatin groups at each time point. Higher osteogenic differentiation was observed in the 2 µM group when compared with the control. Mineralized extracellular deposits were visible after Alizarin red-S staining, and higher mineralization was noted in the 2 and 6 µM lovastatin groups when compared with the 0 µM control. The relative mineralization values of the 0, 2 and 6 µM groups on day 14 were 39.0±9.6, 69.3±6.0 and 60.9±7.5, respectively. This study demonstrated that the application of lovastatin enhanced the osteogenic differentiation of cell spheroids formed from GDSCs. This suggests that combinations of lovastatin and stem cell spheroids may have the potential for use in tissue engineering.

Introduction
The majority of cell culture experiments are conducted on two-dimensional surfaces, including micro-well plates, tissue culture flasks and Petri dishes, due to the ease and convenience of two-dimensional cultures (1). However, two-dimensional cultures may have limitations for the evaluation of cell and tissue physiology, including the communication between a cell and its matrix and between adjacent cells (2). To overcome these limitations, three-dimensional culture techniques have been applied (3,4). Three-dimensional cultures have advan- (3,4). Three-dimensional cultures have advan-. Three-dimensional cultures have advantages including the ability to represent in vivo morphologies and the potential for use in drug discovery with primary and stem cells (5). In a previous study, three-dimensional culture platforms were made using engineered microenvironments, such as highly porous biomimetic scaffolds, which exhibited higher cell differentiation efficiency compared with their two-dimensional counterparts (4). Furthermore, three-dimensional cultures have been shown to support the long-term expansion of nephrogenic progenitor cells (6). Gingiva-derived stem cells (GDSCs) display multipotency with high proliferation characteristics (7,8) Lovastatin is a cholesterol-lowering agent (9); it is involved in regulation of the mevalonate pathway and also affects Akt pathways that are involved in cell proliferation and apoptosis, leading to antiproliferative effect (10). However, the effects of lovastatin on mesenchymal stem cells with three-dimensional cultures have not been well elucidated. Therefore, the purpose of the present study was to evaluate the effects of lovastatin on the proliferation and osteogenic differentiation of human gingiva-derived stem cells (GDSCs) using concave microwells. To the best of the authors' knowledge, this investigation is the first to elucidate the effects of lovastatin on three-dimensional spheroid cultures using mesenchymal stem cells derived from gingiva.

Materials and methods
Isolation and culture of human GDSCs. Gingival

Lovastatin increases the proliferation and osteoblastic differentiation of human gingiva-derived stem cells in three-dimensional cultures
Medicine (no. KC11SISI0348). Informed consent was obtained from all participants according to the Act on Legal Codes for Biomedical Ethics and Safety and the Declaration of Helsinki. Human GDSCs were isolated and cultivated following the protocol published in the present authors' previous study (7). The gingival tissues were collected and maintained in sterile phosphate-buffered saline (PBS; Welgene, Inc.) containing 100 U/ml penicillin and 100 µg/ml streptomycin (Sigma-Aldrich; Merck KGaA) at 4˚C. The tissues were de-epithelialized, separated into 1-2-mm 2 fragments, 0.2 µm filtered, and digested in modified in α-minimum essential medium (α-MEM; Gibco; Thermo Fisher Scientific, Inc.) containing dispase (1 mg/ml; Sigma-Aldrich; Merck KGaA) and collagenase type IV (2 mg/ml; Sigma-Aldrich; Merck KGaA) at 37˚C for 30 min. The cell suspension was filtered with a 70-µm cell strainer (Falcon; BD Biosciences) and the cells were then incubated at 37˚C in a humidified incubator with 5% CO 2 . After 24 h, the non-adherent cells were washed with PBS.
Formation of spheres and evaluation of cellular morphology. Fig. 1 demonstrates the overview of the present study design. Cells were plated onto silicon elastomer-based concave microwells (StemFIT 3D; MicroFIT) of 600 µm diameter at a density of 4x10 5 cells/well and cultured in osteogenic media (StemPro ® Osteogenesis Differentiation Kit; Gibco; Thermo Fisher Scientific, Inc.) at 37˚C. The medium was refreshed at 3-day intervals. To examine the effect of lovastatin, the cells were cultured in the presence of lovastatin (Abcam) at final concentrations of 0 (untreated control), 2 and 6 µM using dimethyl sulfoxide (DMSO) as the vehicle at plating. The concentrations of lovastatin used in the present study were based on those used in previously published studies (11)(12)(13)(14). Equal amounts of DMSO were added to each culture sample to offset the influence of this dissolving vehicle (15). The cells expressed CD44 surface marker and the cell spheroids were positive for SSEA-4 (7,16). The morphology of the microspheres was viewed under an inverted microscope (CKX41; Olympus Corporation) on days 2, 7 and 14 following plating. The diameter of the cell spheroids was measured at each time point.
Determination of cytotoxicity. The cytotoxicity of lovastatin was evaluated on days 2, 7 and 14 with a Cell Counting kit-8 (CCK-8; Dojindo Molecular Technologies, Inc.) according to the manufacturer's protocol. The absorbance at 450 nm was measured spectrophotometrically using a microplate reader (BioTek Instruments, Inc.).
Alkaline phosphatase activity assays. Stem cell spheroids grown with osteogenic media (StemPro ® Osteogenesis Differentiation kit; Gibco; Thermo Fisher Scientific, Inc.) were obtained on day 2 and 8. Alkaline phosphatase activity assays were tested using a commercially available kit (K412-500, BioVision, Inc., Milpitas, CA, USA) following the manufacturer's protocol. The cells were resuspended with an assay buffer, sonicated and then centrifuged at 15,000 x g for 10 min at 4˚C to remove insoluble material. The supernatant was mixed with p-nitrophenylphosphate substrate and incubated at 25˚C for 60 min. The optical density was determined spectrophotometrically at 405 nm.
Alizarin red-S staining. To investigate mineralized nodule formation, cells were grown with osteogenic media (StemPro ® Osteogenesis Differentiation kit; Gibco; Thermo Fisher Scientific, Inc.) for 14 days. The cell spheroids were fixed with 4% paraformaldehyde at room temperature and stained with Alizarin red-S (ScienCell Research Laboratories, Inc.) at room temperature for 30 min. Inverted microscopy (CKX41) was used for evaluation of the stained cells (magnification, x100). The relative value of mineralization was determined by measuring the relative intensity of staining using image processing and analysis software (ImageJ version 1.8.0; National Institutes of Health).
Statistical analysis. Data are presented as the mean ± standard deviation with 95% confidence of intervals (95% CI). Experiments were performed at least three times. A test of normality was performed with a Shapiro-Wilk test. Two-way analysis of variance (ANOVA) was performed for evaluation of the effects of concentration and time and one-way ANOVA was used to determine the differences among groups, followed by Tukey's post hoc test. The analysis was conducted with SPSS 12 for Windows (SPSS, Inc.). P<0.05 was considered to indicate a statistically significant result.

Results
Evaluation of cell morphology and cellular viability. GDSCs formed spheroids in the concave microwells. The morphology of the spheroids at day 2 is shown in Fig. 2A-C. The morphologies of the spheroids at days 7 and 14 were similar to those at day 2 ( Fig. 2D-I). No obvious changes in morphology were observed as the incubation time increased. The diameters of the spheroids were smallest in the 6 µM group at day 2 (P<0.001; Fig. 3 Table I). In general, the diameters of the spheroids were maintained throughout the incubation period.
The CCK-8 assay results for the 0, 2 and 6 µM groups on  Table II).
Alkaline phosphatase activity assay. The results of the alkaline phosphatase activity assay on days 2 and 8 are shown in Fig. 5. The absorbance values at 405 nm on day 2 for the   Mineralization assay. Mineralized extracellular deposits were observed after Alizarin red-S staining on day 14 (Fig. 6).
Higher mineralization was observed in the 2 and 6 µM groups

Discussion
The present study clearly demonstrates that lovastatin at the tested concentrations did not adversely affect the viability of the stem cell spheroids, and increased their osteogenic differentiation.
Statins are drugs that are widely used for lowering serum cholesterol, but have also been shown to enhance new bone formation in vitro and in rodents in previous studies (15,(17)(18)(19). The locally delivery of lovastatin using biodegradable polymer nanobeads of poly(lactic-co-glycolide acid) has been found to improve fracture healing in rats (20). In another study, lovastatin-loaded biodegradable polyurethanes exhibited sustained release of biologically active lovastatin, and the released lovastatin significantly enhanced the osteogenic differentiation of osteoblastic cells in vitro (21). Similarly, the present study using cell spheroids without scaffold confirmed the increased osteogenic differentiation of stem cells when cultured with 2-6 µM lovastatin, suggesting a potential application in therapeutic agents for bone formation.
In a previous study, statins impaired the survival of primary human mesenchymal progenitor cells via mevalonate depletion, nuclear factor κB signaling and B-cell lymphoma 2/adenovirus E1B 19 kDa protein-interacting protein 3 when 1 and 10 µM simvastatin or atorvastatin was used (22). In another study, exposure to simvastatin (0-20 µM) induced a reduction in sphere-forming capacity and cell viability, Table II    Values are presented as the mean ± SD (95% confidence interval). accompanied by a concentration-and time-dependent increase in caspase-3/7 activity (23). Conversely, lovastatin (0.01-1 µM) was found to prevent mesenchymal stem cells from undergoing hypoxia/serum deprivation-induced apoptosis through inhibition of the mitochondrial apoptotic pathway, leading to attenuation of caspase-3 activation (24). In an in vitro study, simvastatin inhibited mesenchymal stem cell apoptosis and increased vascular endothelial growth factor, and combined treatment with simvastatin and mesenchymal stem cells induced a significant improvement in blood reperfusion and a notable increase in capillary density (25). Mesenchymal stem cells have been applied in tissue engineering, for tissues including bone, cartilage, fat and other connective tissue (26). Mesenchymal stem cells have been characterized from a variety of dental-related tissues, including periodontal ligaments, papilla, follicle, dental pulp of exfoliated deciduous and adult teeth, and the maxillary sinus membrane, which represent rich sources of mesenchymal stem cells (27,28). These stem cells have the capacity for self-renewal and multi-lineage differentiation, including osteogenic, chondrogenic and adipogenic differentiation (29). In addi- (29). In addi-. In addition, dental stem cells display several advantages, including a high proliferation rate, high viability and easy induction to distinct cell lineages (27). Moreover, human GDSCs can be harvested during routine practice under local anesthesia and may be considered an excellent source for tissue-engineering purposes (7). Further studies are warranted to evaluate the effects of combination therapy using animal models.
The present study demonstrated that cell spheroids formed from stem cells combined with the application of lovastatin at the tested concentrations had enhanced osteogenic differentiation capability. Therefore, combinations of lovastatin and stem cell spheroids may potentially be useful for tissue engineering purposes.

Availability of data and materials
All data generated or analyzed during the present study are included in the published article.

Authors' contributions
BK, JT, YK and JP designed the study, performed the experi-and JP designed the study, performed the experiments, were responsible for data collection and analysis, and participated in drafting the manuscript. All the authors read and approved the final version of the manuscript.

Ethics approval and consent to participate
All procedures involving human participants were in accordance with the 1964 Helsinki Declaration and its later amendments or comparable ethical standards, and were approved by the Institutional Review Board of the Catholic University of Korea, College of Medicine (no. KC11SISI0348). Informed consent was obtained from the participant.

Patient consent for publication
Not applicable.