Science 12.08.2022_downmagaz.net

RESEARCH | RESEARCH ARTICLES cleavage site (DPRRA) (44, 45) (fig. S8A), affinity and breadth through somatic muta- unmask the fusion peptide region and elicit suggesting that these mAbs can bind before tions suggest that their elicitation may re- S1 shedding or S fusogenic refolding to post- quire multiple rounds of selection, possibly by broadly protective antibody responses. Nota- fusion state, consistent with the preferential heterologous coronavirus infections. A complex binding to prefusion S2 relative to postfusion developmental pathway has been reported for bly, the fusion peptide region was also found S2 (fig. S2A). However, the 2P (K986P, V987P) HIV-1–neutralizing mAbs and may be a general prefusion-stabilizing mutations (4, 5) that lie requirement for mAbs that recognize highly to stimulate broadly reactive CD4 T cells (69), outside of the epitope (Fig. 4A and fig. S8A) conserved epitopes (53, 54). providing a cue for intramolecular help in the inhibited ACE2-mediated enhancement of mAb binding, implying an impediment of receptor- Previous studies have identified serum anti- generation of such antibodies. In conclusion, induced allosteric conformational changes, bodies binding to the fusion peptide of SARS- consistent with recent findings (46). Enhanced CoV-2 and showed, through depletion or our study identifies neutralizing mAbs with binding of fusion peptide–specific mAbs was peptide inhibition experiments, that such also observed with SARS-CoV S– and MERS- antibodies can contribute to the serum- unprecedented breadth and sheds light on a CoV S–expressing HEK293T cells in the pres- neutralizing activity in a polyclonal setting ence of their respective receptors ACE2 (47) (34–38). We show here that some fusion receptor-triggered conformational change that and DPP4 (48). By contrast, fusion peptide- peptide–specific mAbs have direct neutraliz- specific mAbs bound efficiently to HEK293T ing activity in vitro against alpha- and beta- limits the immunogenicity of this conserved cells displaying the alphacoronaviruses NL63 coronaviruses and can ameliorate pathology and 229E S independently of receptor engage- and viral burden in vivo at high doses. Al- region, with potential impacts on universal ment by ACE2 (49) or APN (50), respectively though the neutralizing activity of these mAbs (Fig. 4A and fig. S8A). We observed that VN01H1 is low when used alone, it is possible that in coronavirus vaccine design. and C77G12 neutralized authentic SARS-CoV-2 the context of a polyclonal response, they may Omicron BA.1 and BA.2 more potently than synergize with other antibodies that favor the REFERENCES AND NOTES Washington-1, likely because of increased acces- exposure of the fusion peptide region, as sibility of the fusion peptide in Omicron S (46). shown here with S2E12. 1. A. N. Vlasova et al., Clin. Infect. Dis. 74, 446–454 (2022). Furthermore, the neutralization potencies of 2. J. A. Lednicky et al., Nature 600, 133–137 (2021). VN01H1 and C77G12 could be further improved The biochemical events associated with S2′ 3. M. A. Tortorici et al., Cell 185, 2279–2291.e17 (2022). by engineering mAbs with smaller formats (such proteolytic processing for SARS-CoV, SARS- 4. A. C. Walls et al., Cell 181, 281–292.e6 (2020). as scFv), suggesting possible limited accessi- CoV-2, and MERS-CoV that are critical for 5. D. Wrapp et al., Science 367, 1260–1263 (2020). bility of the fusion peptide region (Fig. 4C). unmasking the fusion peptide region have 6. L. Piccoli et al., Cell 183, 1024–1042.e21 (2020). been described previously (55–59). The mAbs 7. A. J. Greaney et al., Cell Host Microbe 29, 44–57.e9 (2021). It was previously shown that certain mAbs isolated in this study provide tools to charac- 8. A. Z. Wec et al., Science 369, 731–736 (2020). can mimic ACE2 binding and trigger fusogenic terize allosteric conformational changes in S 9. D. Pinto et al., Nature 583, 290–295 (2020). activity (51, 52). Indeed, the addition of S2E12, triggered by receptor binding or by receptor- 10. T. N. Starr et al., Nature 597, 97–102 (2021). an ACE2-mimicking mAb, but not S2M11, a mimicking mAbs (51, 52), a process that has 11. C. A. Jette et al., Cell Rep. 36, 109760 (2021). mAb that locks SARS-CoV-2 S trimer in the thus far eluded structural definition. These 12. M. A. Tortorici et al., Nature 597, 103–108 (2021). closed state (51), enhanced binding of the events lead to exposure of the fusion peptide 13. Y.-J. Park et al., Science 375, 449–454 (2022). fusion peptide-specific mAbs (Fig. 4B). Consistent region for proteolytic cleavage by TMPRSS2 14. D. R. Martinez et al., Sci. Transl. Med. 14, eabj7125 (2022). with this finding, we observed a synergy in or endosomal cathepsins preceding membrane 15. D. Pinto et al., Science 373, 1109–1116 (2021). pseudovirus neutralization when C77G12 was fusion (32). The finding that this region is 16. C. Wang et al., Nat. Commun. 12, 1715 (2021). used in combination with S2E12 (Fig. 4D), readily accessible in alphacoronaviruses is 17. G. Song et al., Nat. Commun. 12, 2938 (2021). suggesting that fusion peptide-specific mAbs may consistent with the higher neutralizing activity 18. M. M. Sauer et al., Nat. Struct. Mol. Biol. 28, 478–486 (2021). be more effective in the context of a polyclonal of some fusion peptide–specific mAbs against 19. P. Zhou et al., Sci. Transl. Med. 14, eabi9215 (2022). antibody response against the RBD. Collec- NL63 and 229E and proposes that conforma- 20. P. Zhou et al., bioRxiv 479488 [Preprint] (2022); tively, these results indicate that the fusion tional changes might occur spontaneously as a peptide epitope only becomes accessible upon consequence of molecular breathing. https://doi.org/10.1101/2022.03.04.479488. receptor-induced S conformational changes in 21. D. Corti et al., Science 333, 850–856 (2011). sarbecovirus SARS-CoV-2, SARS-CoV, and merbe- Broadly neutralizing mAbs recognizing 22. M. M. Sajadi et al., Cell 173, 1783–1795.e14 (2018). covirus MERS-CoV, whereas it is more readily the fusion peptide region of influenza virus 23. H.-X. Liao et al., Nature 496, 469–476 (2013). accessible in the alphacoronaviruses NL63 and hemagglutinin and HIV-1 gp120 have guided 24. J. F. Scheid et al., Science 333, 1633–1637 (2011). 229E, possibly because of molecular breathing the design of universal vaccine candidates against 25. J. Huang et al., Nature 491, 406–412 (2012). to different extents in these S trimers. these highly variable pathogens (60–66). Sim- 26. C. Dreyfus et al., Science 337, 1343–1348 (2012). ilarly, the conservation of the fusion peptide 27. D. Sok, D. R. Burton, Nat. Immunol. 19, 1179–1188 (2018). Discussion region across the Orthocoronavirinae sub- 28. D. Corti, A. Lanzavecchia, Annu. Rev. Immunol. 31, 705–742 (2013). family provides a rationale to develop a vac- 29. D. Pinna, D. Corti, D. Jarrossay, F. Sallusto, A. Lanzavecchia, Using a high-throughput screening of memory cine capable of inducing fusion peptide–specific B cells from immune donors with a panel of responses that may be less prone to viral escape Eur. J. Immunol. 39, 1260–1270 (2009). S proteins, we identified a new class of mAbs because of the potential fitness loss of escape 30. A. C. Walls et al., Cell 185, 872–880.e3 (2022). that target the fusion peptide region. These mutants (67, 68). The finding that prefusion 31. E. M. Anderson et al., Cell 184, 1858–1864.e10 (2021). mAbs have different binding modalities, some stabilizing 2P mutations (4, 5) abolish receptor- 32. M. Hoffmann et al., Cell 181, 271–280.e8 (2020). of which exhibit unprecedented cross-reactivity induced exposure of the fusion peptide suggests 33. C.-L. Hsieh et al., Science 369, 1501–1505 (2020). to the fusion peptide of all four coronavirus that current mRNA vaccines may be suboptimal 34. E. Shrock et al., Science 370, eabd4250 (2020). genera. The findings that the UCAs of some of in inducing antibody responses against this 35. K. W. Ng et al., Science 370, 1339–1343 (2020). these mAbs bind preferentially to common region. The broadly neutralizing coronavirus 36. C. M. Poh et al., Nat. Commun. 11, 2806 (2020). cold coronaviruses and that the mAbs acquire mAbs described in this study may be used as 37. C. Daniel et al., J. Virol. 67, 1185–1194 (1993). probes to design immunogens that can better 38. H. Zhang et al., J. Virol. 78, 6938–6945 (2004). 39. A. C. Walls et al., Nature 531, 114–117 (2016). 40. Q. Xiong et al., bioRxiv 477490 [Preprint] (2022); https://doi.org/10.1101/2022.01.24.477490. 41. F. A. Lempp et al., Nature 598, 342–347 (2021). 42. H. V. Dang et al., Nat. Struct. Mol. Biol. 26, 980–987 (2019). 43. A. C. Walls et al., Proc. Natl. Acad. Sci. U.S.A. 114, 11157–11162 (2017). 44. B. A. Johnson et al., Nature 591, 293–299 (2021). 45. M. Hoffmann, H. Kleine-Weber, S. Pöhlmann, Mol. Cell 78, 779–784.e5 (2020). 46. S. M.-C. Gobeil et al., Mol. Cell 82, 2050–2068.e6 (2022). 47. W. Li et al., Nature 426, 450–454 (2003). 48. V. S. Raj et al., Nature 495, 251–254 (2013). 49. H. Hofmann et al., Proc. Natl. Acad. Sci. U.S.A. 102, 7988–7993 (2005). 50. C. L. Yeager et al., Nature 357, 420–422 (1992). 51. M. A. Tortorici et al., Science 370, 950–957 (2020). 52. A. C. Walls et al., Cell 176, 1026–1039.e15 (2019). 53. D. R. Burton, Nat. Rev. Immunol. 19, 77–78 (2019). 54. M. Bonsignori et al., Immunol. Rev. 275, 145–160 (2017). 55. J.-E. Park et al., Proc. Natl. Acad. Sci. U.S.A. 113, 12262–12267 (2016). 56. J. K. Millet, G. R. Whittaker, Proc. Natl. Acad. Sci. U.S.A. 111, 15214–15219 (2014). 57. P. V. Raghuvamsi et al., eLife 10, e63646 (2021). 58. S. Yu et al., Proc. Natl. Acad. Sci. U.S.A. 119, e2111199119 (2022). 59. E. Qing et al., Cell Rep. 39, 110786 (2022). 60. C. Cheng et al., PLOS ONE 14, e0215163 (2019). 61. K. Xu et al., Nat. Med. 24, 857–867 (2018). SCIENCE science.org 12 AUGUST 2022 • VOL 377 ISSUE 6607 741

RESEARCH | RESEARCH ARTICLES 62. R. Kong et al., Science 352, 828–833 (2016). PHYSICAL CHEMISTRY 63. L. Ou et al., Sci. Rep. 10, 3032 (2020). 64. R. Nachbagauer et al., Nat. Med. 27, 106–114 (2021). Proton-coupled energy transfer in molecular triads 65. H. M. Yassine et al., Nat. Med. 21, 1065–1070 (2015). 66. A. Impagliazzo et al., Science 349, 1301–1306 (2015). Belinda Pettersson Rimgard1†, Zhen Tao2†‡, Giovanny A. Parada2,3, Laura F. Cotter2, 67. S. Belouzard, V. C. Chu, G. R. Whittaker, Proc. Natl. Acad. Sharon Hammes-Schiffer2*, James M. Mayer2*, Leif Hammarström1* Sci. U.S.A. 106, 5871–5876 (2009). We experimentally discovered and theoretically analyzed a photochemical mechanism, which we 68. I. G. Madu, S. L. Roth, S. Belouzard, G. R. Whittaker, J. Virol. term proton-coupled energy transfer (PCEnT). A series of anthracene-phenol-pyridine triads formed a local excited anthracene state after light excitation at a wavelength of ~400 nanometers (nm), 83, 7411–7421 (2009). which led to fluorescence around 550 nm from the phenol-pyridine unit. Direct excitation of phenol- 69. J. S. Low et al., Science 372, 1336–1341 (2021). pyridine would have required ~330-nm light, but the coupled proton transfer within the phenol- 70. M. N. Prichard, C. Shipman Jr., Antiviral Res. 14, 181–205 (1990). pyridine unit lowered its excited-state energy so that it could accept excitation energy from anthracene. Singlet-singlet energy transfer thus occurred despite the lack of spectral overlap ACKNOWLEDGMENTS between the anthracene fluorescence and the phenol-pyridine absorption. Moreover, theoretical calculations indicated negligible charge transfer between the anthracene and phenol-pyridine units. We thank all study participants who donated blood and devoted We construe PCEnT as an elementary reaction of possible relevance to biological systems and time to our research; M. Biggiogero, A. Franzetti Pellanda, future photonic devices. E. Picciocchi, T. Terrot, S. Tettamanti, T. Urbani, L. Vicari, and all personnel at the hospitals and nursing homes for providing blood C oupling of proton transfer to electronic In all the above reactions, the coupled trans- samples; D. Vaqueirinho, S. Jovic, I. Giacchetto Sasselli, transitions exerts a great influence on fer of the proton changes the energy landscape R. Schmidt, and X. Xi from the Sallusto laboratory for help the energetics and kinetics of these of the overall process, which obviously affects with blood processing; M. Kopf (ETH Zurich) for providing processes. Photoacids are well-studied its dynamics and rate. There are also important the pLVX-puro-ACE2 transfer plasmid; H. Tani (University examples, where light excitation leads effects of the much greater mass of the proton, of Toyama) for providing the reagents necessary for preparing VSV to a rapid proton release and, in some cases, whose wave function is more localized than pseudotyped viruses; M. Weisshaar (ETH Zurich) for his artistic even produces long-lived pH changes and pH that of the electrons involved, which may rendition of fig. S1; and Y. Z. Tan from the Rubinstein laboratory for gradients (1, 2). Another example is excited- impose a greater sense of directionality on the helpful discussion. Funding: The work was supported by the state intramolecular proton transfer (ESIPT), process. Understanding and analyzing such Henry Krenter Foundation (F.S.); the European Research Council which generates an excited tautomer state of effects, theoretically and experimentally, is a (AdG no. 885539 ENGRAB to A.L.); the National Institute of Allergy the light-absorbing molecule; systems where fundamentally interesting challenge and of and Infectious Diseases of the National Institutes of Health this state is fluorescent have applications for practical importance. (grants DP1AI158186 and HHSN272201700059C to D.V.); optoelectronic materials (3). Often, ESIPT a Pew Biomedical Scholars Award (D.V.); Investigators in the instead allows for ultrafast dissipation of In 2019, our groups used a series of anthracene- Pathogenesis of Infectious Disease Awards from the Burroughs excess electronic energy as heat, and is there- phenol-pyridine (An-PhOH-Py) triads (Fig. 1) (16) Wellcome Fund (D.V.); Fast Grants (D.V.); the Natural fore an important photostabilizing mecha- to demonstrate the first example of concerted Sciences and Engineering Research Council of Canada (M.M.); nism in plastics (4) as well as in biomolecules PCET occurring in the Marcus inverted region University of Washington Arnold and Mabel Beckman cryoEM such as DNA (5). (17). This behavior, where the rate decreases center; the National Institutes of Health (grant S10OD032290 with increasing driving force, is important for to D.V.); Beamline 5.0.1 at the Advanced Light Source at More recently, the charge transfer mecha- slowing recombination after photochemical Lawrence Berkley National Laboratory; and EOC research nism of proton-coupled electron transfer (PCET) charge separation, but was previously consid- funds. D.V. is an investigator of the Howard Hughes Medical has generated great interest (6). PCET is a key ered unlikely for concerted PCET reactions Institute. F.S. and the Institute for Research in Biomedicine are reaction behind biological energy conversion because of the nearly barrierless tunneling to supported by the Helmut Horten Foundation. Competing in photosynthesis, respiration, and nitrogen higher vibrational product states (18). None- interests: J.S.L., J.J., F.S., A.L., and A.Ca. are currently listed as fixation (7), as well as in synthetic chemical theless, photoexcitation of the An subunit of inventors on multiple patent applications, which disclose the systems for photoredox catalysis and renewable the triads results in the formation of a local subject matter described in this manuscript. A.L., D.C., F.S., J.N., energy conversion and storage (8–11). It is also excited state (LES) on the anthracene (Eq. 1) M.M.-R., L.A.P., G.S., and D.P. hold shares in Vir Biotechnology. The involved in DNA synthesis and repair, photo- (16, 19). In tenths of picoseconds, the LES Veesler laboratory and the Sallusto laboratory have received receptor function, and many other radical en- forms the charge-separated state (CSS) by elec- sponsored research agreements from Vir Biotechnology Inc. The zyme reactions (12–15). Of particular interest tron transfer from the PhOH group to the 1*An remaining authors declare no competing interests. Author is the concerted PCET mechanism, where the concerted with proton transfer to the pyridine contributions: Conceptualization: J.S.L., M.A.T., A.L., D.V., F.S.; electron and proton are transferred without Funding acquisition: F.S., D.V., A.L.; Investigation: J.S.L., J.J., the formation of high-energy intermediates, M.A.T., M.M., D.P., A.C., S.B., J.E.B., A.J., C.S., M.R., A.C.W., D.M., thereby offering a low-barrier reaction pathway. F.M., P.P., D.J., M.F., R.A., B.W., J.No., J.Ne., M.M., L.A.P., C.G., P.F., A.Ce.; Methodology: J.S.L., J.J., M.A.T., M.M., D.P., A.C., 1Department of Chemistry - Ångström Laboratory, Uppsala Fig. 1. Structures of the studied anthracene- S.B., J.E.B., A.J., A.C.W., D.M., F.M., P.P., D.J., M.F., R.A., B.W., University, SE 75120 Uppsala, Sweden. 2Department of phenol-pyridine ([An-PhOH-py]) triads. The enol J.No.; Supervision: J.S.L., D.C., D.V., A.L., F.S.; Writing, review and Chemistry, Yale University, New Haven, CT 06520, USA. form is the stable ground-state form, and the editing: J.S.L., M.A.T., G.S., D.C., D.V., A.L., F.S. Data and materials 3Department of Chemistry, College of New Jersey, Ewing, NJ electronically excited keto form is the local electron- availability: All data associated with this manuscript are available in 08628, USA. proton transfer (LEPT) state. the main text or the supplementary materials, including the FACS data gating strategy. The crystal structures were deposited to the †These authors contributed equally to this work. ‡Present address: Protein Data Bank with accession numbers listed in table S3. All Department of Chemistry, University of Pennsylvania, Philadelphia, further relevant source data that support the findings of this study are available from the corresponding authors upon reasonable request. PA 19104, USA. Materials are available through materials transfer agreements (UBMTAs or similar agreements) with the Institute for Research in *Corresponding author. Email: [email protected] Biomedicine and the University of Washington. License information: This work is licensed under a Creative Commons Attribution 4.0 (S.H.-S.); [email protected] (J.M.M.); International (CC BY 4.0) license, which permits unrestricted use, distribution, and reproduction in any medium, provided the [email protected] (L.H.) original work is properly cited. To view a copy of this license, visit https://creativecommons.org/licenses/by/4.0/. This license does not apply to figures/photos/artwork or other content included in the article that is credited to a third party; obtain authorization from the rights holder before using such material. SUPPLEMENTARY MATERIALS science.org/doi/10.1126/science.abq2679 Materials and Methods Figs. S1 to S8 Tables S1 to S5 References (71Ð106) MDAR Reproducibility Checklist Submitted 28 March 2022; accepted 3 July 2022 10.1126/science.abq2679 742 12 AUGUST 2022 ¥ VOL 377 ISSUE 6607 science.org SCIENCE

RESEARCH | RESEARCH ARTICLES (Eq. 2). The charges subsequently recombine to has transferred to the pyridine. The LEPT negligible charge transfer between the An and the ground-state reactants (Eq. 3), once again state is the same excited tautomeric state as PhOH-py units, which is in sharp contrast to in a concerted electron-proton transfer reac- that formed by ESIPT after excitation of the PCET reactions. This result significantly ex- tion but with a large driving force (>2.0 eV in separate phenol-pyridine compound (22), and pands the scope of concerted processes where CH3CN). This process consequently places the similar to other ESIPT compounds (3). This proton transfer controls reactivity. Many excited- reaction in the inverted region, as was shown state was proposed to be populated after ex- state charge and energy transfer reactions in by the rate dependence on driving force and citation in 4 to 8 (20, 21); it can possibly be proteins and DNA occur between hydrogen- solvent polarity. populated from the CSS (Eq. 4) or be formed bonded pigments, making PCEnT possibly im- directly from the LES (Eq. 5) by a process portant in natural systems. Yet such a reaction An‐PhOH‐py →hn 1ÃAn‐PhOH‐pyðLESÞ ð1Þ here denoted proton-coupled energy trans- has, to the best of our knowledge, not pre- fer, PCEnT (see below). viously been reported. 1ÃAn‐PhOH‐py→An•À‐PhO•‐pyHþðCSSÞ ð2Þ An•À‐PhO•‐pyHþ→An‐ ýPhO‐pyHŠðLEPTÞ Spectroscopic observation of LEPT state An•À‐PhO•‐pyHþ→An‐PhOH‐py ð3Þ ð4Þ The lack of observed CSS formation in triads The CSS intermediate was observed only in 4 to 8, together with the computationally triads 1 to 3, however; in triads 4 to 8, a 1ÃAn‐PhOH‐py P→CEnT An‐ýPhO‐pyHŠðLEPTÞ predicted excited-state keto tautomer (LEPT) similar quenching of the LES did not result in (20), laid the basis for further investigation. an observable CSS (16). This result is puzzling, ð5Þ The triad absorption spectra agree well with a as the charge recombination driving force is sum of contributions from the separate com- larger for 4 to 8, meaning the reaction would Here, we present experimental and theoretical ponents: a vibrationally structured anthracene be more deeply in the inverted region and evidence for direct formation of the fluores- band around 400 nm and a PhOH-py band thus even slower than for 1 to 3. Theoretical cent LEPT state from the LES at 77 K in a around 330 nm (Fig. 2B and fig. S1). The LEPT investigations have suggested the involve- butyronitrile glass, in which formation of the state should fluoresce around 500 nm in ment of another state that would lie lower CSS is inhibited. In this experimental example aprotic, rigid media, as has been shown for the than the CSS in 4 to 8 (20, 21). This is the so- of PCEnT, excitation energy was transferred excited keto state, *[PhO-pyH] of the separate called local electron-proton transfer (LEPT) from 1*An to 1*[PhO-pyH], coupled to proton phenol-pyridine compound (22). Therefore, we state, where the PhOH-py fragment is in a transfer within the PhOH-py unit (Eq. 5). This examined the triads (1, 2, 4 to 7) by fluores- local, electronically excited tautomeric (keto) finding was surprising, as it is a singlet-singlet cence spectroscopy in rigid butyronitrile glass state, *[PhO-pyH], in which the phenolic proton transition without any detectable spectral over- at T = 77 K. We also studied these triads com- lap between the donor fluorescence and accep- putationally and found that the minimum- tor absorption. Moreover, PCEnT occurred with Fig. 2. Steady-state A B absorption (room temperature) and fluores- cence (77 K) spectra in butyronitrile. (A) Sche- matic Jablonski diagram of the possible transitions at 77 K for the triads: radiative decay (downward solid line) and nonradiative decay. The CSS energy depends strongly on the chemical environment (16, 19) and is higher C D than the LES at 77 K. (B) Absorption (solid line) and fluorescence spectra normalized to the first An emission peak (dashed line, lexc: 400 nm) of 1 (black) and the reference compounds of similar concentrations (4 to 7 mM): 2,4-di-tert-butyl-6-(pyridin- 2-yl)phenol (PhOH-py, blue) and 9-cyano-10-meth- ylanthracene (CN-MeAn, green). (C) Fluorescence spectra of triads 1, 2, 4, 5, 6, and 7 with excitation at 400 nm (1 and 2), 365 nm (6), and 375 nm (5 and 7). Note that 7 was measured with twice the absorption at the excitation wavelength. (D) Absorption spectrum (black line), fluorescence excitation spectrum (lem = 455 and 550 nm, solid line) and fluorescence spectrum (lexc = 400 nm, dashed line) of 2. SCIENCE science.org 12 AUGUST 2022 • VOL 377 ISSUE 6607 743

RESEARCH | RESEARCH ARTICLES energy LEPT geometry exhibits a ~90° twist latter is expected for the LEPT state (22) and triads. The relative intensity of the 1*An and between the phenol and pyridine fragments, was similarly observed for all measured triads LEPT fluorescence was independent of triad whereas the planar LEPT geometry is a first- concentration (fig. S2), and control experi- order saddle point. This twisting, which is with an electron-donating group on the pyridine: ments with 9-CN-10-Me-An showed only the known to quench the fluorescence (22, 23), 1, 2, 5, and 6 (Fig. 2, B and C). The LEPT emission structured 1*An fluorescence (Fig. 2B). This is expected to be hindered in the glass at 77 K. shifted to lower energy when the electron- excluded any excimer or other anthracene Thus, the computational studies were per- donating group was weaker (2 > 1 > 6), con- complex as responsible for the broad fluores- formed for both the planar and twisted LEPT sistent with the computational results in cence around 550 nm. geometries. tables S8 and S12. Because 1 and 5 share the same functional group (R2 = 4-Me) on the Fluorescence excitation spectroscopy pro- The 77 K fluorescence spectrum of 1, after pyridine, the LEPT energy should be similar. vided unambiguous evidence for LEPT state selective excitation of the An unit at 400 nm, Still, the peak for 5 was blueshifted by 11 nm, formation from excited 1*An (Fig. 2D). The showed a structured 1*An fluorescence at 420 which indicated a small but significant cou- shape of the excitation spectrum of the 550-nm to 580 nm, accompanied by a broad band with LEPT emission peak (purple line) matched that a maximum around 550 nm (Fig. 2B). The pling to the An unit, which differs in its func- of the absorption spectrum (black line), with tional group (CN in 1; H in 5) between these AB C Fig. 3. Transient fluorescence emission spectra in butyronitrile at 77 K. ~30 ps (instrumental response function IRFFWHM = 65 ps) and decreases (A) Contour map of the transient fluorescence emission intensity counts (dark simultaneously with the the rise of LEPT fluorescence that maximizes after blue = highest intensity) of 6 in butyronitrile (0.2 mM) between 390 and ~140 ps. (C) Time-resolved fluorescence traces (lines, fit; circles, data) within 640 nm within a 2-ns time window. (B) Emission spectra at given times after the LES (445 nm, blue line) and LEPT (570 nm, green line) bands; the IRF is excitation with 370-nm light. The 1*An fluorescence reaches a maximum at shown as a dashed line. AB Fig. 4. Potential energy surfaces and optimized structures for the electronic glass, although a partial twist could occur to lower the energy of the LEPT state states discussed. (A) Schematic potential energy surfaces for the ESIPT and relative to the LES (table S11). GS, ground state. (B) The optimized structures and PCEnT reactions, color-coded to mark the electronically excited fragment the natural bond orbital [NBO (29)] charges on the anthracene (blue), phenol (green, An; purple, PhOH-py). Note that ESIPT is barrierless on the excited-state (red), and pyridine (purple) subunits of 1 for the GS, LES, CSS, and planar (not potential energy surface (purple), whereas PCEnT follows the reorganization twisted, first-order saddle point) LEPT state in the gas phase. The charges reflect of the heavy nuclei including the solvent (green to purple). In fluid solution, the the electron density, where a more negative number equals a more electron-dense LEPT state will twist along the PhOH-py bond, where 90° leads to a conical fragment (the methylene spacer between anthracene and phenol is not included in any intersection, and nonradiative decay. This twist is impeded in the 77 K butyronitrile of the fragments, and therefore those charges are excluded in this analysis). 744 12 AUGUST 2022 ¥ VOL 377 ISSUE 6607 science.org SCIENCE

RESEARCH | RESEARCH ARTICLES contributions from both the PhOH-py (330 nm) fast transient fluorescence spectroscopy mea- The photophysics of triads 1 to 8 have been and the An (300 to 430 nm) absorption bands. surements of 6 with excitation at 370 nm, in previously studied by transient absorption (TA) The excitation spectrum of the 450-nm 1*An the mid- to red part of the An absorption band at ambient temperature (16), and 6 was also emission (green line) showed instead a weaker (Fig. 3). The observed features matched those examined by time-resolved fluorescence spec- relative contribution from the 330-nm PhOH-py seen in the steady-state experiments, and global troscopy at various temperatures (25). In the band. The 550-nm LEPT state could thus be fits were performed within the peak regions TA study at room temperature, the initial an- formed either by exciting the PhOH-py unit, of the LES (380 to 490 nm) and LEPT (520 to thracene LESs in 1 to 3 decayed to a long-lived which would cause formation of the LEPT 655 nm) fluorescence (Fig. 3B). The rise time state (up to t = 755 ps) with spectral features state by the well-established ESIPT mecha- of the LEPT peak was fitted to ~44 ps and a in agreement with the CSS being formed by nism, or directly via excitation of the An unit similarly short decay component could be PCET, including obvious features of the an- followed by energy transfer from 1*An. found for the LES (46 ps), which reinforced the thracene anion. The same mechanism was notion that the LEPT state could be populated assigned to the decay of the LES of 4 to 8 on Evidence for proton-coupled energy directly from the LES (as depicted in the overall the basis of the close correlation of LES decay transfer mechanism reaction scheme in Fig. 2A). rates with the free energies of PCET for each triad. These triads did not show a long-lived The excitation energy transfer from the an- We can exclude LEPT formation via a sequen- transient in the TA spectra at ambient tem- thracene donor to the phenol-pyridine unit tial LES-CSS-LEPT reaction at 77 K, because peratures. However, computations later placed to yield the LEPT state is a type of mecha- the CSS would be substantially destabilized in the LEPT state below the CSS for 4 and 6 (20), nism that we denote proton-coupled energy a solvent glass (Fig. 2A) as a result of its greater indicating that rapid CSS-to-LEPT conversion transfer, PCEnT (Fig. 2A). Energy transfer di- dipole moment relative to both the LES and and/or a contribution of PCEnT from the LES rectly to the *[PhOH-py] state would be strongly LEPT (Fig. 4B). Moreover, the solvent glass at would be reasonable. uphill, as is seen from its low-wavelength ab- 77 K prevented the large degree of twisting of sorption spectrum around 330 nm, but proton the An group to reach the CSS optimized The previously reported spectrally resolved transfer lowers the phenol-pyridine excited-state structure (Fig. 4B and table S1), which further ultrafast fluorescence spectra of 6 at ambient energy to make the PCEnT process energetically destabilized the CSS energy. The CSS lies less temperatures showed an H/D kinetic isotope feasible. The quantum yield for the LEPT for- than 0.2 eV below the LES at room temper- effect (KIE) = 4 and did not show any evidence mation was deemed to be similar irrespective ature (16), and small-molecule charge separa- of emission from a LEPT state. At lower tem- of which fragment was excited, on the basis of tion is typically destabilized by more than peratures, in a 2-methyl-tetrahydrofuran the good agreement between the shapes of the 0.5 eV by freezing out solvent repolarization (2-Me-THF) glass or a polymethylmethacrylate absorption spectrum (Fig. 2D) and the 550-nm in polar solvents (24). Our gas-phase calcula- (PMMA) film, the subnanosecond integrated excitation spectrum, where the respective tions in table S6 indeed placed the CSS higher fluorescence followed a single exponential with PhOH-py and An amplitudes were fairly equal than the LES. Moreover, the LEPT fluores- t ≈ 150 ps, with an unexplained long–time scale in the two spectra. This fact shows that PCEnT cence was formed faster in 6 (~44 ps) than in (>> ns) component. The data in 2-Me-THF from 1*An was about as efficient as ESIPT from 1 (~69 ps), despite greater room-temperature were obtained in the presence of 3% methanol, *[PhOH-py]. driving force for CSS formation (by ~0.5 eV) in which was observed to quench the LEPT emis- the latter (16); this is clear experimental proof sion in the current study. Therefore, the two To follow the formation and decay of the against a reaction involving the CSS. studies are not directly comparable. In retrospect, LEPT state in real time, we performed ultra- Fig. 5. Avoided-crossing region between the S1 and S2 adiabatic electronic virtual orbitals that characterize the S1 excited state, are shown for these states for triad 6, corresponding to a crossing between the LES and LEPT three configurations. For (A), the NTOs are localized on the anthracene (LES), whereas diabatic electronic states, and analysis of the S1 excited state for relevant for (C), the NTOs are localized on the phenol-pyridine fragment (LEPT state), configurations. The energy curves were obtained from an excited-state adiabatic where the proton has mostly transferred. For (B), an approximately equal mixture molecular dynamics trajectory and are plotted along a unitless general reaction of LES and LEPT is observed, as indicated by the two NTOs. Above the orbital coordinate rather than time; this trajectory is far from equilibrium and is not directly drawings of each configuration, the numbers are the charges (in a.u., the elementary comparable to experiment because of the reaction conditions (supplementary text). charge of an electron) of the anthracene, phenol, and pyridine fragments for the S1 (A to C) Configurations obtained from this trajectory. The dominant natural state, in that order (excluding the charge of the CH2 bridging group and the transition orbitals (NTOs), corresponding to the transitions from occupied to transferring hydrogen, which has a charge of ~0.5 a.u.). SCIENCE science.org 12 AUGUST 2022 • VOL 377 ISSUE 6607 745

RESEARCH | RESEARCH ARTICLES these data could be consistent with the LES different proton vibrational states remain to could also form the basis of a molecular switch decaying either by PCET or by the new PCEnT be further explored. pathway reported here in butyronitrile glasses; for energy transfer, where proton transfer additional studies are needed. Figure 5 demonstrates the PCEnT process observed in a molecular dynamics trajectory could be controlled independently by hydrogen- Physical principles underlying proton-coupled for triad 6 propagated on the S1 state, which energy transfer changed from LES to LEPT character as the bonding or (de)protonation, thus effectively proton transferred from the phenol to the pyr- The observation that PCEnT directly formed idine. Representative configurations along this turning energy transfer on or off. We speculate the excited tautomeric LEPT state was sur- trajectory, along with an analysis of the frag- prising and needs mechanistic clarification. As ment charges and dominant natural transition that these processes could already be oper- both the 1*An and the LEPT state are fluores- orbitals, show that the character of the S1 state cent, PCEnT is clearly a singlet-singlet process. changed from LES to a mixture of LES and ational in natural systems with hydrogen- However, there was no detectable spectral LEPT and then to LEPT. The nonadiabatic overlap between the donor fluorescence and coupling between the S1 and S2 adiabatic bonded light-absorbers, such as porphyrins in acceptor ground-state absorption (Fig. 2B), electronic states (27, 28) exhibited a peak at photosynthetic assemblies, DNA and flavin which is a requirement for Förster (dipole- the avoided crossing, providing evidence that proteins, but have yet to be identified. dipole) energy transfer. The absence of low- this process is nonadiabatic. Throughout this energy absorption of the keto form [PhO–-pyH+] process, the charge on the anthracene fragment REFERENCES AND NOTES is because the ground state, at the LEPT ge- remained virtually zero, indicating that the ometry, lies ~0.5 eV higher than the enolic LES-to-LEPT transition did not involve charge 1. Y. Liao, Acc. Chem. Res. 50, 1956–1964 (2017). ground state (table S4) and was thus not transfer between anthracene and the phenol- 2. L. Schulte, W. White, L. A. Renna, S. Ardo, Joule 5, 2380–2394 populated to any appreciable extent. Instead, pyridine group, consistent with the PCEnT as the reaction occurred at quite a small mechanism. Moreover, the net charge trans- (2021). donor-acceptor distance (~6 Å between the fer from the phenol to the pyridine was also 3. A. P. Demchenko, K.-C. Tang, P.-T. Chou, Chem. Soc. Rev. 42, central rings of An and PhOH), the coupl- nearly zero after accounting for the charge of ing was presumably dominated by exchange the transferring proton. In the framework of 1379–1408 (2013). and/or penetration terms (26). Any stron- the vibronically nonadiabatic PCET theory 4. J. Catalan et al., J. Am. Chem. Soc. 114, 5039–5048 ger coupling mechanism was probably pre- (6), the transferring proton would be treated vented by the methylene spacer between the quantum mechanically and would tunnel (1992). subunits. during the LES-to-LEPT transition. 5. A. L. Sobolewski, W. Domcke, C. Hättig, Proc. Natl. Acad. Sci. U.S.A. We propose, therefore, that PCEnT could be Conclusions and outlook 102, 17903–17906 (2005). described by a nonadiabatic surface crossing 6. S. Hammes-Schiffer, A. A. Stuchebrukhov, Chem. Rev. 110, (Fig. 4A), in analogy to Marcus-type theoretical We have experimentally proved the involve- descriptions of PCET (6, 7). Thermal fluctua- ment of the excited *[PhO-pyH] (LEPT) state, 6939–6960 (2010). tions of heavy nuclei lead to a configuration in the photochemical reactions of these triad 7. A. Migliore, N. F. Polizzi, M. J. Therien, D. N. Beratan, where electronic energy transfer and proton molecules. This mechanism competes with tunneling occur with energy conservation PCET from the LES and could explain why Chem. Rev. 114, 3381–3465 (2014). (crossing point of the reactant and product the CSS was not observed in 4 to 8. The LEPT 8. J. J. Warren, T. A. Tronic, J. M. Mayer, Chem. Rev. 110, free energy curves). It is different from PCET, state was formed by a novel proton-coupled however, in that there is no charge transfer energy transfer (PCEnT) mechanism. This is 6961–7001 (2010). between donor and acceptor. The product a singlet-singlet exchange energy transfer, 9. P. R. D. Murray et al., Chem. Rev. 122, 2017–2291 state is a fluorescent excited state when phenol- coupled to proton transfer in the acceptor pyridine twisting is hindered, and can thus be unit, that can be described as a nonadiabatic (2022). readily detected. Moreover, this mechanism transition to the excited product state. Its free 10. N. S. Lewis, D. G. Nocera, Proc. Natl. Acad. Sci. U.S.A. 103, differs from an ESIPT reaction, where proton energy barrier in turn is determined by the transfer is driven by intramolecular charge reorganization energy associated with the 15729–15735 (2006). redistribution, and the electronic coupling is solvent and other heavy nuclei, similar to 11. R. Tyburski, T. Liu, S. D. Glover, L. Hammarström, J. Am. usually strong, which leads to ultrafast reac- Marcus ET theory or PCET theories (6, 17). In tions. The PCEnT reaction described here was contrast to PCET, there is negligible charge Chem. Soc. 143, 560–576 (2021). instead driven by a remote, weakly coupled transfer between donor and acceptor. 12. J. Stubbe, D. G. Nocera, C. S. Yee, M. C. Y. Chang, Chem. Rev. electronically excited state. This mechanism can have several intriguing 103, 2167–2201 (2003). The rearrangement from the enol to the keto implications for both natural and synthetic 13. C. Aubert, M. H. Vos, P. Mathis, A. P. M. Eker, K. Brettel, form of the phenol-pyridine acceptor unit is systems. PCEnT may allow for the photosen- much greater than just the transfer of a proton, sitization of otherwise dark states, such as Nature 405, 586–590 (2000). as it involves substantial charge redistribution seen here by the LEPT formation from the low 14. F. Lacombat et al., J. Am. Chem. Soc. 141, 13394–13409 and changes in bond lengths and angles. These energy excitation of the LES, rather than the combined effects change the excitation energy direct excitation of the PhOH-py unit at higher (2019). of the phenol-pyridine acceptor unit and allow energy (ca 330 nm). Moreover, as energy trans- 15. S. Hu, A. V. Soudackov, S. Hammes-Schiffer, J. P. Klinman, it to come into resonance with that of the fer occurs over long distances but proton anthracene. This behavior is the fundamental transfer is short-range, we suggest that this ACS Catal. 7, 3569–3574 (2017). basis for coupling of proton transfer to energy mechanism may be of importance to direct 16. G. A. Parada et al., Science 364, 471–475 (2019). transfer and distinguishes PCEnT from typi- energy migration in e.g., light-harvesting sys- 17. R. Marcus, N. Sutin, Biochim. Biophys. Acta Rev. Bioenerg. 811, cal energy transfer reactions. The details of tems and photonic materials, where the PCEnT this coupling and the role of coupling with acceptor would form a low-energy site. PCEnT 265–322 (1985). 18. S. J. Edwards, A. V. Soudackov, S. Hammes-Schiffer, J. Phys. Chem. B 113, 14545–14548 (2009). 19. L. F. Cotter, B. Pettersson Rimgard, G. A. Parada, J. M. Mayer, L. Hammarström, J. Phys. Chem. A 125, 7670–7684 (2021). 20. E. R. Sayfutyarova, S. Hammes-Schiffer, J. Am. Chem. Soc. 142, 487–494 (2020). 21. E. R. Sayfutyarova, S. Hammes-Schiffer, J. Phys. Chem. Lett. 11, 7109–7115 (2020). 22. D. LeGourriérec, V. Kharlanov, R. G. Brown, W. Rettig, J. Photochem. Photobiol. Chem. 117, 209–216 (1998). 23. S. Kim, J. Seo, S. Y. Park, J. Photochem. Photobiol. Chem. 191, 19–24 (2007). 24. G. L. GainesIII, M. P. O’Neil, W. A. Svec, M. P. Niemczyk, M. R. Wasielewski, J. Am. Chem. Soc. 113, 719–721 (1991). 25. M. A. Bowring et al., J. Am. Chem. Soc. 140, 7449–7452 (2018). 26. G. D. Scholes, K. P. Ghiggino, J. Phys. Chem. 98, 4580–4590 (1994). 27. S. Fatehi, E. Alguire, Y. Shao, J. E. Subotnik, J. Chem. Phys. 135, 234105–234121 (2011). 28. Q. Ou, G. D. Bellchambers, F. Furche, J. E. Subotnik, J. Chem. Phys. 142, 064114–14 (2015). 29. A. E. Reed, L. A. Curtiss, F. Weinhold, Chem. Rev. 88, 899–926 (1988). 30. B. Pettersson Rimgard et al., Data for “Proton-coupled energy transfer in molecular triads,” Open Science Framework. DOI: 10.17605/OSF.IO/JY46X. ACKNOWLEDGMENTS We acknowledge N. Kaul, L. D’Amario, and E. Sayfutyarova for their expertise and helpful discussions. Funding: Supported by Swedish Research Council grant 2020-05246 and NIH grants 2R01GM50422 and R35GM139449. Author contributions: B.P.R. and L.H. conceived the project and, together with Z.T. and S.H.-S., they wrote the paper. B.P.R., L.F.C., and G.A.P. collected the data, and B.P.R. analyzed all the data. Z.T. and S.H.S. performed 746 12 AUGUST 2022 • VOL 377 ISSUE 6607 science.org SCIENCE

RESEARCH | RESEARCH ARTICLES all the computational work. All authors participated in active the authors, some rights reserved; exclusive licensee American Supplementary Text discussions and reviewed the manuscript. Competing interests: Association for the Advancement of Science. No claim to original Figs. S1 to S11 The authors have no competing interests. Data and materials US government works. www.science.org/about/science-licenses- Tables S1 to S13 availability: All data needed to evaluate the conclusions in the journal-article-reuse References (31Ð50) paper are present in the paper or the supplementary materials. An Data S1 Excel file with steady-state and time-resolved absorption and SUPPLEMENTARY MATERIALS fluorescence data, as well as the optimized geometry for triads science.org/doi/10.1126/science.abq5173 Submitted 13 April 2022; accepted 7 July 2022 1 and 6, and their NBO charges, has been deposited in the Open Materials and Methods Published online 21 July 2022 Science Framework (30). License information: Copyright © 2022 10.1126/science.abq5173 PLANT SCIENCE Agrobacterium from expressing the reporter and confounding fluorescence measurements Synthetic genetic circuits as a means of (fig. S1). Nine of the 10 synthetic activators reprogramming plant roots were able to turn on their target promoters in N. benthamiana (factor of 3 to 45 increase Jennifer A. N. Brophy1,2*, Katie J. Magallon1, Lina Duan1, Vivian Zhong2, Prashanth Ramachandran1, relative to promoter-only controls), and seven Kiril Kniazev1, José R. Dinneny1* of these were specific to the promoters they were designed to regulate (Fig. 1, C and D). The shape of a plantÕs root system influences its ability to reach essential nutrients in the Specificity is needed to construct circuits free soil and to acquire water during drought. Progress in engineering plant roots to optimize of cross-talk between components (20). In con- water and nutrient acquisition has been limited by our capacity to design and build genetic trast to the activators, only four of the synthetic programs that alter root growth in a predictable manner. We developed a collection of synthetic repressors generated more than a factor of 2 transcriptional regulators for plants that can be compiled to create genetic circuits. These circuits change in gene expression, and none achieved control gene expression by performing Boolean logic operations and can be used to predictably alter complete repression (Fig. 1, E and F). Thus, root structure. This work demonstrates the potential of synthetic genetic circuits to control gene many of the repressors required further opti- expression across tissues and reprogram plant growth. mization before they could be used to construct circuits. R ecent advances in our understanding of Constructing synthetic genetic circuits the molecular mechanisms driving plant in plants To achieve greater control over expression, development have identified key regu- To construct genetic circuits in plants, we first we modified the synthetic TFs. Replacing the lators of agronomic traits, but our limited generated a collection of synthetic transcrip- human herpesvirus–derived VP16 activation ability to control gene expression in tional regulators. We designed simple tran- domain with the activation domain from plants is a barrier to applying this knowledge scriptional activators [composed of bacterial Arabidopsis ETHYLENE RESPONSE FACTOR (1, 2). Indeed, it remains challenging to express DNA binding proteins, VP16-activation do- 2 (ERF2AD) improved activity of the PhlF-based genes in specific patterns in plants, especially mains, and SV40 nuclear localization signals synthetic activator by nearly a factor of 6 (Fig. if those patterns are outside the reach of (NLSs)], and synthetic repressors that rely on 1G) (21). A small fragment of ERF2 (78 amino native promoters (3–6). One attractive solu- steric hindrance to achieve repression (com- acids) was used to reduce the possibility of tion is synthetic genetic circuits, which have posed of only DNA binding proteins and NLSs) regulation by native plant components, such been constructed in a wide array of prokary- (Fig. 1, A and B) (12–16). Similar to previous as small RNAs, kinases, etc. (21). In contrast to otic organisms and eukaryotic cell lines (7–10). synthetic promoter designs, activatable plant the activators, adding a repressor domain to However, circuit technology has been difficult promoters were created by fusing six copies of the AmtR-based synthetic repressor either to implement in plants because of the long the DNA sequence (operator) bound by these had no effect [TOPLESS (TPLRD)] or only time scales required for producing transgenic transcription factors (TFs) to a minimal plant modestly increased repression (relative to lines and the difficulty of tuning circuit activ- promoter [positions –66 to +18 of the cauli- NLS only) by 20% [ETHYLENE RESPONSE ity across heterogeneous cell types (11). Here, flower mosaic virus (CaMV) 35S promoter] FACTOR 4 (ERF4RD)] or 33% [SALT TOL- we show that quantitative transient expres- (Fig. 1A) (17, 18). Repressible promoters were ERANCE ZINC FINGER (ZAT10RD)] (Fig. sion assays can be used to aid the construction built by placing one operator sequence at the 1H) (22, 23); thus, swapping the activation of circuits in plants by providing a platform for 3′ end of a full-length CaMV 35S promoter domain led to increased expression, whereas testing and tuning circuit designs. We show (Fig. 1B) (19). This design was selected to avoid adding repressive domains had a more modest how these assays can be leveraged to improve disrupting 35S promoter activity when adding effect. the performance of circuits designed to pro- operators. duce specific spatial patterns of gene expres- We also engineered the synthetic promoters sion across root tissues, and we use one of Synthetic TF activity was measured in plants to tune expression levels and increase dynamic the circuits to precisely modify plant root by delivering DNA encoding the proteins and range. Changing the number of TF binding sites architecture. their target promoters to Nicotiana benthamiana in the AmtR-activatable promoter resulted in leaves using Agrobacterium tumefaciens. To a collection of promoters of varying strength. 1Department of Biology, Stanford University, Stanford, CA, enable quantitative measurements of gene ex- Promoters with three or more operators am- USA. 2Department of Bioengineering, Stanford University, pression, we used the synthetic promoters to plified expression relative to the constitutive Stanford, CA, USA. drive expression of green fluorescent protein G1090 promoter used to drive expression of *Corresponding author. Email: [email protected] (J.A.N.B.); (GFP) and normalized GFP expression to a con- AmtR-VP16 (Fig. 1I, compare 6× to 3× with [email protected] (J.R.D.) stitutively expressed mCherry encoded on the control) (24). In contrast, promoters with two same T-DNA (see supplementary materials). operators roughly replicated G1090’s expres- We also introduced an intron to GFP to prevent sion level and one operator produced less expression (Fig. 1I). These data suggest that synthetic promoters can be used to tune gene expression levels by either amplifying or SCIENCE science.org 12 AUGUST 2022 • VOL 377 ISSUE 6607 747

RESEARCH | RESEARCH ARTICLES Fig. 1. Basic building blocks for constructing synthetic genetic circuits. containing one (1×) to six (6×) AmtR operators co-infiltrated with a plasmid (A and B) Schematics of the synthetic transcriptional activators (A) and encoding proG0190::AmtR-VP16. Control construct (C) is proG0190::GFP. repressors (B) built to control gene expression in plants. Small bent arrow (J) Engineered synthetic promoters’ response to AmtR-NLS. In all panels, expression denotes the transcription start site. Activity (C) and specificity (D) of the is measured in N. benthamiana leaves. TFs are expressed using the constitutive synthetic activators. Yellow box highlights orthogonal activators. (E) Activity G1090 promoter. Box plots show the median of 12 leaf punches collected from three of the synthetic repressors. (F) Fluorescence of leaves infiltrated with a leaves that were infiltrated and measured on different days. Box plot hinges plasmid that does not encode GFP. (G) Comparison of PhlF-based activators indicate the first and third quartiles. Dots show individual data points. Stars containing either the VP16 (V16) or ERF2 (E2) activation domain. (H) Comparison denote significant differences in normalized GFP expression and letters denote of AmtR-based repressors containing either no repressive domain or a significance groups (P < 0.01, Student’s two-tailed t test); n.s., not significant. repressive domain from TOPLESS (TPL), ZAT10 (Z10), or ERF4 (E4). Relative change for significantly active transcription factors was calculated by Control (−) is promoter alone. (I) Activity of synthetic activatable promoters dividing the average ON state by the average OFF state. dampening the transcriptional signal produced and PhlF DNA binding proteins, which dem- to include repressive domains, such as ERF4RD, by an input promoter. onstrated strong activation and repression in to prevent synthetic activators from initiat- our initial TF designs, served as the inputs to ing transcription at composite promoters In line with a previous study, we found that all circuits, and GFP served as the output. Cir- (1-1 state, A NIMPLY B; Fig. 2C and fig. S3). moving the location of the operator to a region cuit activity was measured in N. benthamiana However, this requirement could be over- between the TATA and CAAT boxes within leaves, which were infiltrated with multiple come by increasing the number and location the 35S promoter improved repression and Agrobacterium strains, each containing one of repressor operators within the synthetic dynamic range (Fig. 1J) (10). Adding a second plasmid that encoded either an input TF or promoter (B NIMPLY A; fig. S3). Most circuits operator between the TATA box and the tran- the output (Fig. 2, A and B). involved several design-build-test cycles, which scription start site further improved repression were facilitated by the rapid N. benthamiana– and dynamic range (Fig. 1J). However, adding Computation was performed by synthetic based assays and the modular nature of the a third binding site at the 3′ end of the promoters that responded to the input TFs in synthetic biology parts generated here. Addi- promoter reduced its ON state without signif- unique ways. Simple gates, like the A and B tional tuning information is available in figs. icantly reducing the repressed state, resulting BUFFER gates, used synthetic promoters that S3 and S4. in a smaller dynamic range than the two- responded to one input TF; more complex gates operator promoter (Fig. 1J, 58× versus 64×). required synthetic promoters that responded to Layered logic gates in which AmtR and These design features were transferrable to multiple inputs (Fig. 2C). To create functional PhlF do not directly control expression of GFP, other synthetic repressor-promoter pairs (fig. circuits, both the promoter architecture and but instead modulate the expression of other S2), demonstrating that the number and loca- the synthetic TFs needed to be optimized synthetic TFs, occasionally encountered prob- tion of operators can be adjusted to tune pro- (fig. S3). We found that the arrangement of lems where expression cassettes were insuffi- moter activity. operators in the OR promoter affected fold ciently insulated from each other. The A IMPLY change, with the best design containing B gate, which worked well when its output Using these synthetic regulators, we con- alternating pairs of operators (Fig. 2C and fig. genes were encoded on separate plasmids, had structed circuits that perform Boolean logic S3). Additionally, synthetic repressors needed an erroneously reduced “no input” state when operations. Synthetic TFs built with the AmtR 748 12 AUGUST 2022 • VOL 377 ISSUE 6607 science.org SCIENCE

RESEARCH | RESEARCH ARTICLES Fig. 2. Logic gates in N. benthamiana leaves. (A) Schematic for testing Fig. 3. Patterning gene expression using logic gates. (A) Cell types in the circuit activity in N. benthamiana using a transient expression assay (see Arabidopsis root tip: columella (C), lateral root cap (R), epidermis (E), cortex (X), supplementary materials). N. benthamiana leaves were infiltrated with endodermis (N), stele (S), quiescent center (Q). (B) Confocal images show the Agrobacterium strains containing plasmids that encode either one input expression pattern of input promoters SOMBRERO (SMB) (input A) and PIN-FORMED transcription factor (A or B) or the output gene GFP. mCherry was encoded 4 (PIN4) (input B). (C) Confocal images (left) and expected expression patterns on the output plasmid and used to normalize fluorescence measurements. (right) of logic gates that use the SMB and PIN4 promoters to express input TFs. (B) Example images of leaves infiltrated with TRUE and FALSE gate Output of each circuit is nuclear localized GUS-GFP fusion protein. In all panels, components. Circles encompass areas of infiltration. (C) Gate behavior in T-DNA schematics show the identity and arrangement of circuit components. Root N. benthamiana. Bar charts show output of each circuit, reported as the images were taken 5 days after sowing. Box-and-whisker plots show quantified ratio of GFP to mCherry. Green bars indicate gate states that should be ON; reporter expression for three T2 plants from a single transgenic line. Green boxes gray bars indicate states that should be OFF to implement correct logic. indicate cell layers that should be ON and gray bars indicate states that should Data are mean and SD of 12 leaf punches collected from three leaves be OFF to implement correct logic. Expression levels of individual cells in each infiltrated and measured on different days. Dots show individual data points. root layer are shown as black dots. Box plot hinges indicate the first and third Data for individual GFP and mCherry channels are provided in the quartiles. Letters denote significant differences in expression (P < 0.01, StudentÕs supplementary materials. two-tailed t test). Additional independent lines in fig S6. both output genes were encoded on the same A IMPLY B output genes (25). However, in- Engineering spatial patterns of plasmid (Fig. 2C and fig. S4). In an attempt to sulation only improved gate activity by a factor gene expression fix this problem, we added one or two copies of 2 (fig. S4, version 2 versus versions 3 and 4). Functional plant circuits were transferred to of the previously characterized insulator Consequently, the other IMPLY gate (B IMPLY the model plant A. thaliana to test their capac- Arabidopsis thaliana MATRIX ATTACH- A) was built with two output plasmids that ity to generate specific spatial patterns of gene MENT REGION 10 (MAR10) between the two were co-delivered to plant cells (Fig. 2C). expression across root tissues. The tissue-specific SCIENCE science.org 12 AUGUST 2022 ¥ VOL 377 ISSUE 6607 749

RESEARCH | RESEARCH ARTICLES Fig. 4. Engineering root branch density using BUFFER gates. (A) Root length. (C) Primary root length measured from root-hypocotyl junction to the structure (top) and root hairs (bottom) of Arabidopsis plants engineered to root tip. Control samples are labeled as follows: slr-1 (slr), proGATA23::slr (G), modify root branch density. BUFFER gate output promoters contain either one wild type (WT). In (B) and (C), box plots show the median of at least 20 T1 (1×) or two (2×) copies of the consensus AmtR operator or one copy of a plants. Box plot hinges indicate the first and third quartiles. Letters denote mutated AmtR operator (M1 to M3). Scale bar, 1 cm. (B) Lateral root density significant differences in LR density or primary root length (P < 0.01, StudentÕs calculated as the number of emerged lateral roots divided by total primary root two-tailed t test). All measurements and images were taken 10 days after sowing. promoters of SOMBRERO (proSMB, expres- Fig. 3C, and fig. S5). This change was made to ses have been difficult to test directly because sed in the entire root cap) and PIN-FORMED eliminate potential interactions between the genetic changes that affect root structure often 4 (proPIN4, expressed in columella, root cap, promoter (proPIN4) and the downstream gene. have pleiotropic effects. For example, a gain-of- and stele) were used to drive expression of For A NIMPLY B, we removed the ERF4RD function mutation in the developmental reg- our input TFs (Fig. 3, A and B) (26, 27). By from the input TF (PhlF-ERF4). ERF4 is ulator INDOLE-3-ACETIC ACID INDUCIBLE combining the activity of these input pro- believed to recruit histone-modifying enzymes 14 (IAA14) called solitary root (slr-1) eliminates moters using logic gates, we expected to generate and may maintain a repressed state in cells root branching, but also hinders root gravit- several different spatial patterns of gene ex- that are not actively expressing the repressor ropism, root hair development, and primary pression (Fig. 3C, cartoons). (23). Successful A NIMPLY B gates were built root growth (Fig. 4 and fig. S8, slr-1) (30). To with repressors that lack the repression do- disentangle root branching from other develop- We tested eight different logic gates in main (PhlF-NLS) (Fig. 3C and fig. S5). Tuning mental processes, we expressed the slr-1 mutant Arabidopsis, four of which generated the ex- of the remaining gates (OR and NOR) and a gene using a tissue-specific promoter that is pected expression pattern upon first attempt description of our troubleshooting process only ON in lateral root stem cells, proGATA2323. (TRUE, FALSE, A BUFFER, NOT A; Fig. 3C). are outlined in figs. S5 to S7. When expressed from proGATA23, slr-1’s impact Successful gates qualitatively matched the on development is restricted to lateral roots. In expected expression patterns and produced After tuning, all gates qualitatively matched these plants, no lateral roots form, but gravit- a significant difference between GFP expres- the expected expression patterns (Fig. 3C). How- ropism, root hair development, and primary sion in tissues expected to be ON versus those ever, three had expression in a single cell layer root growth are normal (Fig. 4 and fig. S8, expected to be OFF. By this definition of that was either aberrantly high (B BUFFER, proGATA23::slr-1). success, TRUE, FALSE, A BUFFER, and NOT endodermis) or low (OR, stele; NOR, QC). Thus, A gates were successful. For NOT A, the quantitative analysis highlights the challenge We designed BUFFER gates to express difference between the lowest ON (cortex) of implementing circuits across cell types in slr-1 at varying levels in lateral root stem cells and highest OFF (root cap) states was only heterogeneous tissue and additional optimi- to determine whether lateral root branch de- 1.2×; this finding suggests that further opti- zation would be required for these gates to velopment could be quantitatively controlled. mization may be necessary for applying the achieve significant differences across every These BUFFER gates use the GATA23 promoter circuit in other contexts (Fig. 3C). tissue layer. to drive expression of the AmtR-VP16 synthetic TF, which then activates a synthetic promoter The gates that did not work revealed differ- Modifying root structure with one, two, four, or six copies of the AmtR ences in circuit behavior between Arabidopsis operator to drive expression of slr-1. Because and N. benthamiana. For example, the first B To demonstrate how precise spatial control plant transformation randomly inserts trans- BUFFER gate produced a spatially expanded over gene expression may be used to engineer genes into the genome, a single construct can expression pattern relative to the input promoter, development, we modified a key aspect of create a range of expression. Although this with aberrant expression in the quiescent center Arabidopsis root architecture: lateral root variation potentially alleviates the need for (QC) and neighboring initial cells (figs. S5 and branch density. Lateral roots allow plants to promoters of varying strength, none of the S6). Similarly, the A NIMPLY B pattern was radially sample soil, and the number of lateral plants (>20 independent lines per construct) incorrect in Arabidopsis; expression was missing roots that a plant generates affects its ability containing proGATA23::slr-1 or a BUFFER in several root cap cell layers (figs. S5 and S6). to search for water and essential nutrients in gate with two, four, and six copies of the AmtR the environment (28). The close relationship operator produced lateral roots; this suggests We tuned the B BUFFER, A NIMPLY B, OR, between root growth and plant fitness has that the range of expression conferred by these and NOR gates to improve performance. For B led to proposals of ideal root architectures constructs was above the threshold for fully BUFFER gate, we changed the input TF from for plant growth in specific environments blocking lateral root development (Fig. 4, A AmtR-VP16 to QacR-VP16, another highly active (called root ideotypes) (29). Ideotype hypothe- synthetic activator in N. benthamiana (Fig. 1C, 750 12 AUGUST 2022 • VOL 377 ISSUE 6607 science.org SCIENCE

RESEARCH and B; data not shown for 4× and 6× pro- 9. X. J. Gao, L. S. Chong, M. S. Kim, M. B. Elowitz, Science 361, 34. A. Khakhar, D. F. Voytas, Front. Plant Sci. 12, 668580 (2021). moters). To further reduce the strength of the 1252–1258 (2018). 35. J. W. Wang et al., Mol. Plant 12, 1037–1040 (2019). output promoter, we mutated key residues in 36. B. Cole et al., Commun. Biol. 4, 962 (2021). the AmtR binding site (fig. S9) (31). Weaker 10. K. A. Schaumberg et al., Nat. Methods 13, 94–100 (2016). 37. S. Vickovic et al., Nat. Commun. 13, 795 (2022). synthetic promoters resulted in BUFFER gates 11. T. K. Kassaw, A. J. Donayre-Torres, M. S. Antunes, K. J. Morey, that produced a gradient of lateral root den- ACKNOWLEDGMENTS sities, roughly correlated with the strength of J. I. Medford, Plant Sci. 273, 13–22 (2018). the synthetic promoter (Fig. 4 and fig. S10). 12. A. J. Meyer, T. H. Segall-Shapiro, E. Glassey, J. Zhang, We thank T. Beekman for the Arabidopsis slr mutant and T. Vellosillo None of the BUFFER gates affected primary for assistance in establishing imaging protocols. Funding: root growth, root hair development, or gravit- C. A. Voigt, Nat. Chem. Biol. 15, 196–204 (2019). Burroughs Wellcome Fund CASI (J.A.N.B.); Chan Zuckerberg Biohub ropism (Fig. 4, B and C, and figs. S8 and S10). 13. B. C. Stanton et al., Nat. Chem. Biol. 10, 99–105 (2014). (J.A.N.B., J.R.D.); US Department of Energy Biological and Thus, BUFFER gates enabled specific tuning 14. I. Sadowski, J. Ma, S. Triezenberg, M. Ptashne, Nature 335, Environmental Research Program grant DE-SC0008769 (J.R.D.); of lateral root density. Our results show that Faculty Scholar grant from Howard Hughes Medical Institute and the difference in regulatory capacity conferred 563–564 (1988). the Simons Foundation 55108515 (J.R.D.). Author contributions: by our synthetic system has a greater effect on 15. B. C. Stanton et al., ACS Synth. Biol. 3, 880–891 (2014). Conceptualization, J.A.N.B., J.R.D.; data acquisition, J.A.N.B., root branching than the T-DNA insertion site, 16. C. Dingwall, R. A. Laskey, Trends Biochem. Sci. 16, 478–481 (1991). K.J.M., L.D., V.Z., K.K.; visualization, J.A.N.B.; funding acquisition, but they do not define how a specific dosage of 17. R. X. Fang, F. Nagy, S. Sivasubramaniam, N. H. Chua, Plant Cell J.A.N.B., J.R.D.; project administration, J.R.D.; supervision, J.A.N.B., slr-1 determines lateral root development; the J.R.D.; writing–original draft, J.A.N.B., J.R.D.; writing–review and resources generated here can be used to ex- 1, 141–150 (1989). editing, J.A.N.B., V.Z., J.R.D. Competing interests: The authors plore this relationship further. 18. L. Laplaze et al., J. Exp. Bot. 56, 2433–2442 (2005). declare that they have no competing interests. Data and 19. J. T. Odell, F. Nagy, N.-H. Chua, Nature 313, 810–812 (1985). materials availability: All plasmid materials and bacterial strains Our work presents a roadmap for the design 20. J. B. Lucks, L. Qi, W. R. Whitaker, A. P. Arkin, Curr. Opin. will be made available through AddGene. Sequence files are and implementation of synthetic genetic cir- available at doi.org/10.5281/zenodo.6949657 and raw data are cuits that program gene expression across cell Microbiol. 11, 567–573 (2008). available as supplementary materials. License information: types in plants. This approach expands the 21. J. Li et al., Plant Biotechnol. J. 11, 671–680 (2013). Copyright © 2022 the authors, some rights reserved; exclusive impact that a handful of characterized promot- 22. E. Pierre-Jerome, S. S. Jang, K. A. Havens, J. L. Nemhauser, licensee American Association for the Advancement of Science. No ers can have on our ability to express genes claim to original US government works. www.science.org/about/ in specific cells and presents an alternative to E. Klavins, Proc. Natl. Acad. Sci. U.S.A. 111, 9407–9412 (2014). science-licenses-journal-article-reuse searching genomes for tissue-specific promoters 23. M. Ohta, K. Matsui, K. Hiratsu, H. Shinshi, M. Ohme-Takagi, with desired expression patterns. The ability to SUPPLEMENTARY MATERIALS modulate expression in a tissue-specific manner Plant Cell 13, 1959–1968 (2001). provides an opportunity to probe gene dosage 24. F. Ishige, M. Takaichi, R. Foster, N.-H. Chua, K. Oeda, Plant J. science.org/doi/10.1126/science.abo4326 effects in a tissue-specific manner, which previ- Materials and Methods ously lacked a reliable framework in plants. As 18, 443–448 (1999). Figs. S1 to S10 circuit technology matures in plants, it could 25. A. Pérez-González, E. Caro, Sci. Rep. 9, 8474 (2019). References (38–54) be applied to reprogram responses to the en- 26. M. Kamiya et al., Development 143, 4063–4072 (2016). Data S1 and S2 vironment or to better understand gene ex- 27. M. D. M. Marquès-Bueno et al., Plant J. 85, 320–333 (2016). pression control in plants, as has been done 28. R. Rellán-Álvarez, G. Lobet, J. R. Dinneny, Annu. Rev. Plant Biol. Submitted 3 February 2022; accepted 23 June 2022 recently in Drosophila (32, 33). 10.1126/science.abo4326 67, 619–642 (2016). Reprogramming crops using synthetic genetic 29. J. P. Lynch, New Phytol. 223, 548–564 (2019). circuits will require careful tuning. As evidenced 30. H. Fukaki, S. Tameda, H. Masuda, M. Tasaka, Plant J. 29, here, controlling expression levels across tis- sues can be challenging, and although some ap- 153–168 (2002). plications may tolerate aberrant expression in a 31. K. Hasselt, S. Rankl, S. Worsch, A. Burkovski, J. Biotechnol. 154, few cell layers, others will be extremely sensi- tive to off-target expression. To efficiently engi- 156–162 (2011). neer genetic circuits that meet specific design 32. J. Crocker, G. R. Ilsley, D. L. Stern, Nat. Genet. 48, 292–298 (2016). requirements, researchers can leverage the 33. J. Crocker, A. Tsai, D. L. Stern, Cell Rep. 18, 287–296 (2017). growing number of transient expression methods available for plants (34, 35). Advances in single- ◥ cell measurements and spatially resolved omics will also support the construction of circuits in REPORTS multicellular organisms by enabling more precise measurements (36, 37). Ultimately, methods for B I O M AT E R I A L S programming novel traits in plants will become increasingly useful as climate challenges grow and Controlled tough bioadhesion mediated by ultrasound new agricultural solutions are needed. Zhenwei Ma1, Claire Bourquard2, Qiman Gao3, Shuaibing Jiang1, Tristan De Iure-Grimmel4, Ran Huo1, REFERENCES AND NOTES Xuan Li1, Zixin He1, Zhen Yang1, Galen Yang5, Yixiang Wang6, Edmond Lam5,7, Zu-hua Gao8, Outi Supponen2*, Jianyu Li1,9* 1. S. J. Gurr, P. J. Rushton, Trends Biotechnol. 23, 275–282 (2005). 2. A. E. Richardson et al., Science 374, 1377–1381 (2021). Tough bioadhesion has important implications in engineering and medicine but remains challenging to 3. S. Ali, W.-C. Kim, Front. Plant Sci. 10, 1433 (2019). form and control. We report an ultrasound (US)Ðmediated strategy to achieve tough bioadhesion 4. J. Bai et al., Plant Biotechnol. J. 18, 668–678 (2020). with controllability and fatigue resistance. Without chemical reaction, the US can amplify the 5. R. Wang et al., Sci. Rep. 5, 18256 (2015). adhesion energy and interfacial fatigue threshold between hydrogels and porcine skin by up to 100 and 6. L. Chen et al., Agriculture 11, 1195 (2021). 10 times. Combined experiments and theoretical modeling suggest that the key mechanism is US- 7. A. A. K. Nielsen et al., Science 352, aac7341 (2016). induced cavitation, which propels and immobilizes anchoring primers into tissues with mitigated barrier 8. D. Mishra et al., Science 373, eaav0780 (2021). effects. Our strategy achieves spatial patterning of tough bioadhesion, on-demand detachment, and transdermal drug delivery. This work expands the material repertoire for tough bioadhesion and enables bioadhesive technologies with high-level controllability. B ioadhesive technologies find frequent particularly critical for bioadhesives based on use in wearable electronics, biomedical physical interactions such as polymer inter- implants, wound management, anasto- penetration; the polymers are too slow and mosis, regenerative medicine, and drug even impossible to diffuse and entangle with delivery (1–5). However, their use has tissues (6), resulting in poor bioadhesion (7). long been hindered by the barrier effects of Chemical strategies have thus far been pri- biological tissues, such as low permeability marily used for tough bioadhesion. Despite and limited functional groups. Skin, for in- achieving high adhesion energy, they do not stance, contains dense stratum corneum, enable high-level control over bioadhesion limiting the penetration and bonding of bio- in space and time. Exceptions require so- adhesive agents (Fig. 1A). These issues are phisticated surface patterning, exogenous SCIENCE science.org 12 AUGUST 2022 ¥ VOL 377 ISSUE 6607 751

RESEARCH | REPORTS chemicals, and an external apparatus to re- trate tissues by passive diffusion (Fig. 1A) (6). replacing it with phosphate-buffered saline move interfacial bonding (8, 9). Other dis- We hypothesized that US could propel and yielded weak bioadhesion, even with US treat- advantages include interference with payloads immobilize the anchoring agents deep into ment (fig. S3). The adhesion performance was for drug delivery (10), low fatigue threshold tissues, which would form a bridging net- dependent on the pH difference between the (caused by limited functional groups on tissue work at the interface when a hydrogel patch is hydrogel patch and the primer solution (fig. S4) surfaces) (11), and acute and/or chronic toxic- applied (Fig. 1B). This strategy of US-mediated but independent of blood exposure on tissue ity (caused by chemical reactions and their bioadhesion differentiates from existing strat- substrates (fig. S5). We also used an ultrasonic reagents) (12). egies based on passive diffusion, chemical scaler (20 to 35 kHz) like the type used in dental reactions, or invasive mechanical interlock- clinics to obtain tough adhesion (~800 J m−2) To form and control tough bioadhesion, ing, including those in medical devices such on porcine skin with US treatment (fig. S6). as suture anchors and swellable microneedles our strategy leverages ultrasound (US) and (17) (fig. S1). Our strategy is applicable to a large rep- ertoire of materials. The adhesion enhance- anchoring agents such as polymers and nano- US-mediated bioadhesion was achieved in ment by US was confirmed with different two steps. First, we used an ultrasonic trans- hydrogels, including a double-network poly(N- particles. US has been extensively used in the ducer (20 kHz) to apply US to a primer solution isopropylacrylamide)-alginate (PNIPAm-alg) clinical setting for imaging (13), monitoring or suspension of anchoring agent spread on a hydrogel and a single-network PAAm hy- (14), tumor ablation (15), and drug delivery tissue substrate (e.g., freshly excised porcine drogel; the variance of adhesion energy can (16). The anchoring agents can form a physi- skin). Second, we covered the treated area be linked with the bulk toughness of the hy- cally cross-linked network but do not pene- with a hydrogel patch, which triggers gelation drogel (fig. S7). We also demonstrated US- of the anchoring agent at the interface (Fig. 1B). mediated bioadhesion with various anchoring 1Department of Mechanical Engineering, McGill University, As a model system, we deployed a chitosan agents, including gelatin, ChsNC, and CHO- Montréal, Quebec H3A 0C3, Canada. 2Institute of Fluid solution and a polyacrylamide-alginate (PAAm- CNC (Fig. 1D), and the different adhesion Dynamics, Department of Mechanical and Process alg) hydrogel as the primer solution and the performances could be related to the charge Engineering, ETH Zürich, 8092 Zürich, Switzerland. 3Faculty hydrogel patch, respectively. We measured the density and structure of anchoring agents. of Dental Medicine and Oral Health Sciences, McGill adhesion energy between the hydrogel and Gelatin could lead to >200 J m−2 of adhesion University, Montréal, Quebec H3A 1G1, Canada. 4Department the tissue with peeling tests (Fig. 1C). With energy with US treatment, markedly higher of Bioengineering, McGill University, Montréal, Quebec H3A US treatment (116 W cm−2, 1 min), the ad- than that without US (~26 J m−2). CHO-CNC 0E9, Canada. 5Department of Chemistry, McGill University, hesion energy obtained on porcine skin was yielded high adhesion energy (~180 J m−2) Montréal, Quebec H3A 0B8, Canada. 6Department of Food ~1750 J m−2, >15 times that of the no-US con- with US treatment, whereas the condition Science and Agricultural Chemistry, McGill University, Sainte- trol (Fig. 1D). The adhesion energy reached without US yielded weak adhesion (~12 J m−2). Anne-De-Bellevue, Quebec H9X 3V9, Canada. 7Aquatic and ~100 J m−2 within a minute and then pla- ChsNC could not penetrate deep into skin Crop Resource Development Research Centre, National teaued in 10 min (fig. S2). The anchoring agent without the aid of US (fig. S8). The adhesion Research Council of Canada, Montréal, Quebec H4P 2R2, is critical for tough bioadhesion because merely of ChsNC without US remained weak even Canada. 8Department of Pathology and Laboratory Medicine, University of British Columbia, Vancouver, BC V6T 1Z7, Canada. 9Department of Biomedical Engineering, McGill University, Montréal, Quebec H3A 2B4, Canada. *Corresponding author. Email: [email protected] (J.L.); [email protected] (O.S.) Fig. 1. Robust and versatile US-mediated tough bioadhesion. (A) Schematic of skin with barrier effects limiting passive diffusion and impairing bioadhesion. (B) US actively propels and anchors primer agents into a tissue substrate, causing spatially confined tough adhesion between hydrogel and tissue. (C) Representative force-displacement curves of hydrogel-tissue (porcine skin) hybrids with or without US treat- ment in peeling tests. (D) US enables diverse anchoring agents for tough bioadhesion on skin. Chi, chitosan; ChsNC, chitosan nanocrystal; CHO-CNC, aldehyde functionalized cellulose nanocrystal. (E) Adhesion enhancement with US on diverse biological tissues, including skin, buccal mucosa, and aorta. Chitosan was used as the anchoring primer. Data are shown as means ± SD for n = 3 independent experiments. 752 12 AUGUST 2022 ¥ VOL 377 ISSUE 6607 science.org SCIENCE

RESEARCH | REPORTS when they formed amide bonds with tissues The tissue-dependent adhesion performance thermal effect, we monitored the temper- can be related to the mechanics and chemistry ature change using a thermal camera. After (fig. S9). When the US was applied, the ad- of specific tissues (18). The versatility of US- a 1-min application of 16 W cm−2 US, the sur- hesion energy increased to 500 J m−2, two mediated bioadhesion and its indispensable face temperature of porcine skin increased role in nanoparticle bioadhesion unlock the by ~10°C (fig. S13). To determine the effects of orders of magnitude higher than the same potential of various materials for tough bio- this, we replicated the temperature change material without US treatment (~5 J m−2). adhesion (fig. S11) (19). with a temperature-controlled oven to incu- bate hydrogel-primer-tissue hybrids. Because In addition to the skin, our strategy is appli- We investigated the mechanism underlying no significant change of bioadhesion was found US-mediated bioadhesion as follows. The US (fig. S14), we excluded the thermal effect. We cable to other biological tissues, including exerts mainly thermal and mechanical effects then focused on the mechanical effects of US on the primer solution and the tissue sub- such as cavitation, viscous stress, and acoustic buccal mucosa and aorta (Fig. 1E). The mea- strate. Given the US intensity used in this radiation force. Of these, cavitation was hypothe- sured adhesion energy was ~295 J m−2 for study, we detected minimal change of the US- sized to play a key role (21) because viscous buccal mucosa and ~297 J m−2 for aorta; the treated chitosan in terms of structure and stress and acoustic radiation force are less gelling behavior (fig. S12) (20). To evaluate the values for no-US conditions were ~12 and ~17 J m−2, respectively. Tough bioadhesion was evidenced by the debonding of buccal membrane from the underlying tissues during peeling (fig. S10), indicative of strong adhesion. Fig. 2. US-induced cavitation regulates bioadhesion. (A) Experimental bubble clouds, interfacial adhesion energy, and shear strength. (I) Schematic setup for characterizing US-induced microbubble cavitation. (B) Digital illustration of the fatigue test of hydrogel-tissue hybrid. Data are shown as images of the microbubble cloud at peak intensity in a cycle captured means ± SD for n = 3 independent experiments. (J) Representative curves of by the high-speed camera. Scale bar, 100 mm. (C) Processed binary images the cycle number and maximum displacement per cycle at two energy release rates of the bubble clouds. (D) Normalized peak intensity of the induced microbubble clouds as a function of the US intensity. (E) Representative for samples with US treatment. (K) Crack extension rate (dc/dN) versus applied force-displacement curves of hydrogel-tissue hybrids in lap-shear tests. energy release rate G = F/W for hydrogel adhesion on tissues after no US or (F and G) Correlation between applied US intensity and interfacial adhesion exposure to low (16 W cm−2)– or high (116 W cm−2)–intensity US. The linear energy (F) and shear strength (G). (H) Correlation between peak intensity of extrapolation to the G-axis (solid lines) gives the fatigue threshold G0. (L) Correlation of US intensity and obtained fatigue threshold Go. SCIENCE science.org 12 AUGUST 2022 ¥ VOL 377 ISSUE 6607 753

RESEARCH | REPORTS C D AB EF p (kPa) Fig. 3. Spatial control of US-mediated bioadhesion. (A) Schematics of US theoretically (blue circle). Data are shown as means ± SD for n = 3 independent application at varying distances d between the US horn and the tissue substrate. experiments. (E) Side view of theoretically computed fields. The area where (B) Correlation of d and the adhesion energy characterized by lap-shear tests. p < pv during part of the acoustic cycle is shown by a red dashed line for d = 1, 3, (C) Images of porcine skin where a clear boundary between tough adhesion and 5 mm. (F) Corresponding top view at substrate level for d = 1, 3, and 5 mm. and nonadhesion regions can be observed. Scale bar, 1 cm. (D) Correlation of d Red arrowhead indicates the boundary of tough adhesion area. Blue areas are and tough adhesion area measured experimentally (green dot) and computed where total pressure p goes below vapor pressure pv according to theory. significant for the low-frequency US used here dynamics with laser-induced cavitation ex- We next showed that the US treatment me- (see the supplementary text and fig. S15). To periments, showing that the laser-induced diates tough bioadhesion. Because the US test this hypothesis, we combined experiments bubble collapses and jets toward the tissue effects scaled with the distance between the and theoretical modeling to study the link be- substrate (fig. S17 and movie S1); US-induced transducer and the tissue (d), simply maneuver- tween cavitation and bioadhesion. microbubble jetting was suggested by pit ing the US transducer could control the bio- generation in aluminum foil (fig. S18). There adhesion in magnitude and space (Fig. 3A). We first characterized cavitation and bio- might also be contributions from other me- The adhesion energy decreased with increasing adhesion as a function of US intensity. The chanical effects of the US such as viscous d and eventually to the level without US when cavitation was manifested with dynamic micro- stress, which requires further investigation. d > 4 mm (Fig. 3B). The distance d also mediated bubble clouds present under the US transducer the area of bioadhesion (Fig. 3C). The tough (Fig. 2, A and B). A high-speed camera revealed To further study the US effect on interfacial adhesion area decreased with the gap be- the pattern and geometry of oscillating vapor bonding, we performed fatigue fracture tests. tween the transducer and the tissue (Fig. 3D). bubble clouds. The peak cavitation intensity Such tests output interfacial fatigue threshold To understand and predict the spatially con- was quantified with the maximum area of or intrinsic work of adhesion G0, which scales trolled bioadhesion, we conducted theoretical microbubble cloud in each cycle (Fig. 2C). The with the density of interfacial bonds according modeling on the acoustic field produced by US substantially enhanced the peak intensity to the Lake-Thomas theory (18). In fatigue tests the US transducer between the horn and the of cavitation, but the effect diminished when under 180° peeling configuration (Fig. 2I), we substrate (see the supplementary text); the the intensity was high (Fig. 2D). A similar varied the magnitudes of cyclic loads applied experimentally measured displacement of trend in interfacial adhesion energy was ob- onto the specimens for different energy release the US transducer was used in the model tained in modified lap-shear tests (Fig. 2E). rates (G) and monitored the crack extension (fig. S19). We extracted the area on the sub- The adhesion energy and shear strength be- over cycles (dc/dN) (Fig. 2J). The crack growth strate where the absolute pressure dropped tween hydrogel and porcine skin increased rate increased with the loading G (Fig. 2K). below the vapor pressure in every acoustic with the US intensity (Fig. 2, F and G), and The linear regression of the G-dc/dN curves cycle, thereby enabling the formation of cav- both were shown to correlate linearly with the informed the intrinsic work of adhesion (G0). itation bubbles (Fig. 3, E and F). At various peak bubble intensity (Fig. 2H). The adhesion The US treatment raised G0 from ~5 J m−2 values of d, we obtained substantially different energy scaled almost linearly with the dura- (no-US) to ~65 J m−2 (Fig. 2L). With the ab- pressure profiles on the tissue substrate from tion of US treatment (fig. S16), further corrob- sence of covalent bonding at the interface, the which the regions affected by cavitation were orating the correlation between cavitation and high G0 exceeded the case of forming inter- estimated (Fig. 3D). The computational results bioadhesion. facial amide bonds through carbodiimide agreed with the experimental measurements. chemistry (~25 J m−2) (11). These results sub- The spatial resolution of bioadhesion could be This correlation can be explained as follows. stantiate the existence of strong interfacial further improved with focused-US technology When the bubbles collapse and jet, the an- interactions, which are typically only observed (22). In addition to demonstrating spatial con- choring agents are propelled and anchored with covalent bonding, in contrast to the often trol, we demonstrated that US could enable into tissue substrates, enabling interfacial weak physical interactions resulting from in- temporal control over bioadhesion by using bonding with the hydrogel patch. In support terdiffusion, such as entanglement. of this point, we visualized the single-bubble 754 12 AUGUST 2022 ¥ VOL 377 ISSUE 6607 science.org SCIENCE

RESEARCH | REPORTS Fig. 4. In vivo assessments. (A) Schematics of the procedures of US-mediated hybrid using chitosan (D) or ChsNC (E) as the anchoring primer, respectively. bioadhesion with a rodent model. (B) Spatial control of bioadhesion on rat skin. Scale bar, 100 mm. (F) Degree of inflammation after a 1-day hydrogel attachment Blue dotted circle indicates the adhesion region. Scale bar, 1 cm. (C) Adhesion mediated by US; 0 indicates no inflammation; 1, very mild inflammation; 2, mild energy between rat skin and hydrogel without or with US treatment (16 W cm−2). inflammation; 3, moderate inflammation; 4, severe inflammation; and 5, very Data are shown as means ± SD for n = 3 independent experiments. (D and severe inflammation. Statistical significance and P values were determined E) Representative hematoxylin and eosin–stained images of the hydrogel-rat skin by two-sided Student’s t test. *P < 0.05; ***P < 0.001. the heating effect of US and a thermo-gelling Here, we report US-mediated bioadhesion 23. C. Bourquard, O. Supponen, Scripts for: Controlled tough gelatin as the anchoring primer (see the sup- to precisely control hydrogel bioadhesion in bioadhesion mediated by ultrasound, Zenodo (2022); plementary text and fig. S20). The controlla- space and time. Along with drug-eluting hy- https://zenodo.org/record/6718010. bility is particularly desired because current drogels, our strategy enables tough bioadhesion bioadhesives are limited by homogeneous and transdermal delivery of protein concur- ACKNOWLEDGMENTS adhesiveness, poor control over the diffusion rently (see the supplementary text and fig. S22). of adhesive agents, and complications asso- The universal applicability of our strategy We thank L. Mongeau, A. Higgins, A. Moores, and D. Kurdyla for ciated with patterning adhesives (2, 8, 9). promises impacts in broad areas ranging from providing materials and equipment; the Life Sciences Complex wearable devices to drug delivery. Advanced BioImaging Facility (ABIF) at McGill University for Finally, we validated the safety and efficacy fluorescence image collection and processing; and X. Yang at the of US-mediated bioadhesion with a rodent REFERENCES AND NOTES University of Kansas for helpful discussions. Funding: This work was model in vivo (Fig. 4A). This animal study is supported by the Natural Sciences and Engineering Research Council critical because the mechanical index for the 1. S. Nam, D. Mooney, Chem. Rev. 121, 11336–11384 (2021). of Canada (grant RGPIN-2018-04146 to J.L. and grant RGPIN-2019- high-frequency US safety criteria for imag- 2. Z. Ma, G. Bao, J. Li, Adv. Mater. 33, e2007663 (2021). 04498 to Y.W.); the Canada Foundation for Innovation (grant 37719 ing applications is not applicable to the low- 3. J. Li et al., Science 357, 378–381 (2017). to JL); Fonds de Recherche du Quebec - Nature et Technologies frequency US used in this study (16). We tested 4. H. Yuk et al., Nature 575, 169–174 (2019). (grant FRQ-NT NC-270740 to J.L. and B2X doctoral research both chitosan and ChsNCs to demonstrate 5. K. Liu et al., Bioact. Mater. 13, 260–268 (2021). scholarship 268520 to Z.M.); the Canada Research Chair Program in vivo applicability. Tough bioadhesion 6. J. Yang, R. Bai, Z. Suo, Adv. Mater. 30, e1800671 (2018). (J.L.); ETH Zürich (O.S. and C.B.); and the National Research Council formed within minutes selectively on a cir- 7. H. R. Brown, T. P. Russell, Macromolecules 29, 798–800 Canada Industrial Biotechnology program (E.L.). Author contributions: cular region treated with US, indicative of Conceptualization: Z.M., J.L.; Investigation: Z.M., C.B., Q.G., S.J., spatially controlled adhesion (Fig. 4B). By (1996). T.D.G., R.H., X.L., Z.H., Z.Y., G.Y., Z.G.; Methodology: Z.M., J.L., C.B., testing freshly excised rat skin, we confirmed 8. A. Mahdavi et al., Proc. Natl. Acad. Sci. U.S.A. 105, 2307–2312 O.S.; Resources: J.L., O.S., Y.W., E.L.; Supervision: J.L., O.S.; Writing – the high adhesion energy achieved with US original draft: Z.M., J.L., C.B., O.S.; Writing – review and editing: Q.G., in vivo (Fig. 4C). Histological assessments (2008). S.J., T.D.G., R.H., X.L., Z.H., Y.W., E.L., Z.G., Z.M., J.L., C.B., O.S. concluded that there was no marked tissue 9. X. Chen, H. Yuk, J. Wu, C. S. Nabzdyk, X. Zhao, Proc. Natl. Competing interests: The authors declare no competing interests. damage or acute inflammation caused by US Data and materials availability: All data are available in the main (16 W cm−2) or bioadhesives (Fig. 4, D to F). Acad. Sci. U.S.A. 117, 15497–15503 (2020). text or the supplementary materials or have been deposited at Zenodo Rodent skin is more sensitive to the heating 10. J. Li, D. J. Mooney, Nat. Rev. Mater. 1, 16071 (2016). (23). License information: Copyright © 2022 the authors, some of high-intensity US compared with porcine 11. X. Ni, C. Chen, J. Li, Extreme Mech. Lett. 34, 100601 (2020). rights reserved; exclusive licensee American Association for the skin (see the supplementary text and fig. S21). 12. N. Artzi, T. Shazly, A. B. Baker, A. Bon, E. R. Edelman, Adv. Advancement of Science. No claim to original US government works. Future study is needed to establish the US https://www.science.org/about/science-licenses-journal-article-reuse conditions (e.g., intensity and profile) for clin- Mater. 21, 3399–3403 (2009). ical translation. 13. C. Demené et al., Nat. Biomed. Eng. 5, 219–228 (2021). SUPPLEMENTARY MATERIALS 14. C. Wang et al., Nat. Biomed. Eng. 2, 687–695 (2018). 15. C. Lovegrove, Nat. Clin. Pract. Oncol. 3, 8–9 (2006). science.org/doi/10.1126/science.abn8699 16. C. M. Schoellhammer et al., Sci. Transl. Med. 7, 310ra168 (2015). Materials and Methods 17. S. Y. Yang et al., Nat. Commun. 4, 1702 (2013). Supplementary Text 18. Z. Yang, Z. Ma, S. Liu, J. Li, Mech. Mater. 157, 103800 (2021). Figs. S1 to S22 19. S. Rose et al., Nature 505, 382–385 (2014). References (24–33) 20. G. Kim et al., Proc. Natl. Acad. Sci. U.S.A. 116, 10214–10222 (2019). Movie S1 21. F. Duck, T. Leighton, J. Acoust. Soc. Am. 144, 2490–2500 (2018). 22. A. R. Rezai et al., Proc. Natl. Acad. Sci. U.S.A. 117, 9180–9182 Submitted 3 January 2022; accepted 22 June 2022 10.1126/science.abn8699 (2020). SCIENCE science.org 12 AUGUST 2022 ¥ VOL 377 ISSUE 6607 755

RESEARCH | REPORTS ORGANIC CHEMISTRY with the exhaustive fluorination of all cage vertexes and the latent instability of (CF)n on Electron in a cube: Synthesis and characterization of account of the overcrowded fluorine atoms perfluorocubane as an electron acceptor (13–15). By contrast, perfluorocubane [(CF)8, 1], whose structure was proposed in 2004 (16), Masafumi Sugiyama1, Midori Akiyama1*, Yuki Yonezawa1, Kenji Komaguchi2, Masahiro Higashi3, should be sufficiently stable for isolation given Kyoko Nozaki1, Takashi Okazoe1,4 that the vicinal fluorine atoms are sterically less hindered than those in larger polyhedranes. Fluorinated analogs of polyhedral hydrocarbons have been predicted to localize an electron within Moreover, a theoretical study by Irikura (8) their cages upon reduction. Here, we report the synthesis and characterization of perfluorocubane, and our own density functional theory (DFT) a stable polyhedral fluorocarbon. The key to the successful synthesis was the efficient introduction of calculations have suggested that 1 can be ex- multiple fluorine atoms to cubane by liquid-phase reaction with fluorine gas. The solid-state structure of pected to exhibit pronounced electron-accepting perfluorocubane was confirmed using x-ray crystallography, and its electron-accepting character was character [electron affinity = 1.6 eV (8); energy corroborated electrochemically and spectroscopically. The radical anion of perfluorocubane was of the lowest unoccupied molecular orbital examined by matrix-isolation electron spin resonance spectroscopy, which revealed that the unpaired (ELUMO) = −2.8 eV; Fig. 1B]. Here, we report electron accepted by perfluorocubane is located predominantly inside the cage. the synthesis, isolation, and characterization of 1 as a polyhedral fluorocarbon that can F or organic chemists, polyhedral mole- become a focal point of attention and the sub- accept and store an electron within its internal cules such as cubane (1), dodecahedrane ject of numerous theoretical studies (Fig. 1A). cubic cavity (i.e., electron in a cube). (2), and buckminsterfullerene (3) repre- One of their most interesting features is their sent attractive synthetic targets beyond electron-accepting character inside their cage, As a synthetic route to 1, the stepwise intro- which arises from a stabilized vacant orbital duction of eight F atoms onto cubane is im- their arguably subjective elegance. After within the cage derived from multiple s* or- practical; in fact, only two F atoms have been bitals of C–F bonds (8–10). This internal local- successfully introduced into a cubane scaffold the syntheses and structural characterizations ization of electrons stands in stark contrast to common p-conjugated electron acceptors, 1Department of Chemistry and Biotechnology, Graduate of these polyhedral cages were achieved, inter- which usually host electrons on their molecu- School of Engineering, The University of Tokyo, Bunkyo-ku, lar surfaces. Although some experimental evi- Tokyo 113-8656, Japan. 2Department of Applied Chemistry, est quickly shifted to the possibility of encap- dence that supports the formation of (CF)n Graduate School of Advanced Science and Engineering, (where n = 20 or 60) has been reported (11, 12), Hiroshima University, Kagamiyama, Higashi-Hiroshima, sulating guests within their internal cavities. isolation as a single isomer has not yet been Hiroshima 739-8527, Japan. 3Department of Molecular achieved because of the difficulties associated Engineering, Graduate School of Engineering, Kyoto Early examples included the encapsulation of University, Nishikyo-ku, Kyoto 615-8510, Japan. 4AGC Inc., a variety of atoms in fullerene (4) and dode- Suehiro-cho, Tsurumi-ku, Yokohama, Kanagawa 230-0045, cahedrane (5) as well as the encapsulation of Japan. a molecule in fullerene (6, 7). Meanwhile, *Corresponding author. Email: [email protected] polyhedral fluorocarbons (CF)n—i.e., the per- fluorinated analogs of polyhedranes—have Fig. 1. A perfluorinated cube-shaped molecule accepts an electron inside its cage. (A) Schematic representation of the encapsulation of electrons, atoms, or molecules inside cage-shaped molecules. Polyhedral fluorocarbons can accept an electron inside their cage, whereas common p-conjugated electron acceptors usually store electrons on their surface. (B) Structure of perfluorocubane (1) and its predicted properties. The energy level of 1 was calculated at the B3LYP/6-311++G(d,p) level of theory, and the predicted electron affinity has been reported previously (8). (C) Previously reported related molecules and synthetic strategy for this study. 756 12 AUGUST 2022 • VOL 377 ISSUE 6607 science.org SCIENCE

RESEARCH | REPORTS Fig. 2. Synthetic routes to fluorinated cubanes. (A) Prepa- ration of heptafluorinated com- pound 5 and octafluorinated 1. The yield for 3 in the direct fluorination represents an average value of six experiments, and the error (±s) was obtained from the standard deviation (s). THF, tetrahydrofuran; LiHMDS, lithium bis(trimethylsilyl)amide; NFSI, N-fluorobenzenesulfonimide; Bn, benzyl. (B) Preparation of hexafluorinated compound 7. so far (Fig. 1C) (17, 18). Several per-substituted Fig. 3. Crystal structures of fluorocubanes. (A) Thermal ellipsoid plots of 1, 4, and 7 (ellipsoid probability, cubanes have been synthesized through the 50%). (B) Packing structure of 1 in the unit cell. Two molecules in the unit cell are depicted using a space-filling model to show the manner of interaction (carbon, white; fluorine, yellow). (C) Hirshfeld surface for dimers [2+2] photocycloaddition of the correspond- of 1 in the crystal structure mapped using dnorm (left), together with a hypothetical illustration of the ing octa-substituted syn-tricyclo[4.2.0.02,5]octa- interactions in the dimers (right). (D and E) Schematic illustration of noncovalent s-hole interactions. 3,7-dienes (19, 20). Although it is feasible to conceive the synthesis of 1 from syn- octafluorotricycloocta-3,7-diene, its conversion to 1 has not yet been reported (21, 22). Against this background, the introduction of multiple fluorine atoms through radical C–H fluorina- tion using fluorine gas might seem the most promising strategy to synthesize 1. However, a reported attempt to fluorinate cubane using fluorine gas resulted in ring opening owing to the highly strained nature of the cubic skele- ton (23). We envisaged that the C–H fluori- nation of cubane could be achieved using a modified version of the liquid-phase direct fluorination designated the “PERFECT” meth- od, in which the direct fluorination of a par- tially fluorinated ester with an excess of fluorine gas in a fluorous solvent affords the corre- sponding perfluorinated product with marked suppression of C–C bond cleavage (24, 25). This synthetic approach enabled the system- atic preparation of octa-, hepta-, and hexa- fluorocubanes, which enabled study of the effect of the degree of fluorination of cubane on its electron affinity. The synthetic route to 1 is shown in Fig. 2A. A solution of cubanemonoester 2 in 1,2-dichloro- 3-(2-chloro-1,1,2,2-tetrafluoroethoxy)-1,1,2,3,3- pentafluoropropane (CFE-419) was treated at −20°C with an excess of fluorine gas. After the reaction, three new peaks were observed in the 19F nuclear magnetic resonance (NMR) spectrum with an area ratio of 3:3:1, which indicated the formation of the heptafluoro- cubane monoester (3) with an average NMR yield of 15% (for details, see the supplementary SCIENCE science.org 12 AUGUST 2022 • VOL 377 ISSUE 6607 757

RESEARCH | REPORTS Fig. 4. Properties of fluorocubanes. (A) UV-vis absorption spectra of 1 (red), 5 (blue), and 7 (gray) in MeCN; the vertical dashed lines represent the theoretically predicted absorption coefficients and wavelengths calculated at the B3LYP/6-311++G(d,p) level of theory. The red asterisk indicates the predicted absorption wavelength of 1, for which the predicted oscillator strength f is 0.00 given that it is symmetry forbidden. (B) LUMOs and the energy levels of frontier orbitals for 1, 5, and 7, calculated at the B3LYP/6-311++G(d,p) level of theory. Red and green areas depict different phases of the plotted orbitals. (C) ESR spectrum of 1 in a hexamethylethane matrix after g-irradiation at 77 K (black) and simulated ESR spectrum (red) using g = 1.9985 and a(19F) = 19.62 mT (8F). Exp, experimental; Sim, simulated. (D) Spin density (contour = 0.005; top) and singly occupied molecular orbital (SOMO) (isovalue = 0.05; bottom) in the radical anion of 1, calculated at the UB3LYP/6-311++G(d,p) level of theory. materials). Ester 3 was subsequently subjected The structures of fluorocubanes 1, 4, and 7 in the Hirshfeld surface analysis (Fig. 3C) (30). to transesterification without isolation, and were unambiguously determined using single- The interaction energy was estimated to be crystal x-ray diffraction analysis (Fig. 3A). The 3.5 kcal/mol using natural energy decompo- benzyl heptafluorocubanecarboxylate (4) was identical lengths (1.570 Å) of the 12 C–C bonds sition analysis (31). This type of interaction successfully isolated. Upon hydrolysis of 4, a in 1 clearly indicate a nondistorted cubic has recently been recognized as a noncovalent subsequent decarboxylation furnished hepta- structure in the solid state, and this length was carbon-bonding interaction, which is a subset almost identical to a previously reported the- of tetrel bonding. Specifically, carbon bonding fluorocubane (5; 76% yield). The decarboxyla- oretical value that was obtained from high- represents an interaction between a Lewis- tion of perfluorinated tertiary carboxylic acids level DFT calculations (27). The bond length acidic s-hole centered on a carbon atom and was also almost identical to that of the parent a Lewis base, akin to halogen bonds (Fig. 3D) is known to occur without requiring any other cubane C8H8 (1.572 Å) (28). According to Bent’s (32). The only reported example of a synthe- rule, atoms tend to direct the bonding of sized carbon-bonding donor is tetracyanocy- reagents (26). The reaction of 5 with lithium hybrid orbitals of greater p character toward clopropane (33). Bauzá et al. have predicted bis(trimethylsilyl)amide (LiHMDS) and N- electronegative substituents, such as fluorine. that 1 can be expected to act as a strong carbon- fluorobenzenesulfonimide (NFSI) gave 1 in As a result, the s character of the C–C bond is bonding donor that interacts with anions such 51% yield. The 19F and the 19F-decoupled 13C increased, and the bond is shortened. By con- as cyanide (Fig. 3E) (34), and the packing trast, the average C–C bond expands with in- structure of 1 in the single crystal obtained NMR spectra of 1 in acetone-d6 each feature creasing numbers of fluorine atoms, which in the present study strongly supports this one singlet peak [19F: −197.19 parts per million is the result of the repulsion between vicinal prediction. (ppm); 13C{19F}: 103.8 ppm], and a cross peak fluorine atoms (29). These two opposing effects appear to cancel each other out, and hence Fluorocubanes 1, 5, and 7 showed ultraviolet- between these two peaks was observed in the 1 shows a C–C bond length that is similar to visible (UV-vis) absorption peaks at wavelengths 19F-13C heteronuclear multiple quantum co- that of the parent cubane (supplementary longer than 160 nm owing to their lower-lying materials). In the crystal structure, 1 exhibits a LUMOs compared with those of other fluoro- herence (HMQC) spectrum. The Fourier trans- characteristic intermolecular interaction. One alkanes, which are transparent at 160 nm (35). fluorine atom of 1 interacts with the center As shown in Fig. 4A, 1 exhibits the longest form infrared (FTIR) spectrum in the range of the cyclobutane ring of the nearest mole- absorption wavelength edge among the three 400 to 4000 cm−1 contained only two peaks cule, as depicted using a space-filling model in fluorocubanes, and the edges of 5 and 7 shift at 1371 and 798 cm−1, which reflects the high Fig. 3B; this interaction is also clearly observed hypsochromically with decreasing numbers of fluorine atoms. DFT calculations at the B3LYP/ symmetry (Oh point group) of 1. In the mi- croscopic Raman spectrum, the characteristic peak corresponding to the breathing vibration mode of 1 was observed at 552 cm−1. Hexa- fluorocubane (7) was also synthesized from cubanediester 6 (32%, in three steps) using a similar procedure (Fig. 2B and the supplemen- tary materials). 758 12 AUGUST 2022 • VOL 377 ISSUE 6607 science.org SCIENCE

RESEARCH | REPORTS 6-311++G(d,p) level of theory revealed that the coupling constant (hfcc), the theoretically 28. R. J. Doedens, P. E. Eaton, E. B. Fleischer, Eur. J. Org. Chem. LUMO of each fluorocubane comprises the s* predicted hyperfine structure characteristic 2017, 2627–2630 (2017). orbitals of the C–F bonds and that the energy of second-order splitting could be confirmed of the LUMO decreases with increasing num- (37). The isotropic line shape was well sim- 29. J. Cioslowski, L. Edgington, B. B. Stefanov, J. Am. Chem. Soc. ulated by assuming identical coupling of the 117, 10381–10384 (1995). bers of fluorine atoms (Fig. 4B). The results of eight equivalent 19F nuclei. This result suggests that the observed radical anion of 1 rotates rap- 30. M. A. Spackman, D. Jayatilaka, CrystEngComm 11, 19–32 time-dependent DFT calculations suggest that idly even in the glass-state matrix at 77 K. The (2009). the absorptions of 5 and 7 should most likely hfcc (19.62 mT) is in accord with the value be attributed to the transition from the doubly based on the calculated structure (19.24 mT), 31. E. D. Glendening, A. Streitwieser, J. Chem. Phys. 100, in which the spin density is mainly distributed 2900–2909 (1994). degenerate the second-highest occupied molec- inside the cubane cage (Fig. 4D). ular orbital (HOMO−1), which consists of the 32. A. Bauzá, T. J. Mooibroek, A. Frontera, Angew. Chem. Int. Ed. C–C s bonds of the cubane scaffold, to the Previous polyhedral molecules have found 52, 12317–12321 (2013). LUMO. The bathochromic shift of the absorp- widespread applications in functional mate- tion peak of 5 (231 nm) compared with that of 7 rials after their synthesis. Thus, our study, 33. A. Frontera, C 6, 60 (2020). (200 nm) is the result of the lower LUMO level which demonstrates that perfluorinated cage 34. A. Bauzá, T. J. Mooibroek, A. Frontera, Phys. Chem. Chem. of 5 relative to that of 7 and consistent with the compounds act as electron acceptors, can be calculated values (5: 247 nm; 7: 211 nm). expected to pave the way for the molecular Phys. 16, 19192–19197 (2014). Perfluorocubane (1) has a triply degenerate design of distinctive functional organic ma- 35. G. Bélanger, P. Sauvageau, C. Sandorfy, Chem. Phys. Lett. 3, HOMO because of its highly symmetric struc- terials [a connection to the environmental accumulation of per- and polyfluoroalkyl sub- 649–651 (1969). ture, and the predicted oscillator strength of stances (PFAS) cannot be discarded at present, 36. A. M. ElSohly, G. S. Tschumper, R. A. Crocombe, and attention should be paid to this issue the HOMO-LUMO transition is 0.00 (Fig. 4B, when devising applications of such fluorinated J. T. Wang, F. Williams, J. Am. Chem. Soc. 127, 10573–10583 cubanes] (38). (2005). red asterisk). This can be rationalized in terms 37. R. W. Fessenden, J. Chem. Phys. 37, 747–750 (1962). of the point group of 1 (Oh), which renders the REFERENCES AND NOTES 38. M. G. Evich et al., Science 375, eabg9065 (2022). transition from the HOMO (t2u) to the LUMO (a1g) symmetry forbidden (16). In the experi- 1. P. E. Eaton, T. W. Cole, J. Am. Chem. Soc. 86, 3157–3158 ACKNOWLEDGMENTS mental UV-vis absorption spectrum of 1, a (1964). weak absorption was observed at a wave- We thank R. Takita (The University of Tokyo) for atmospheric 2. R. J. Ternansky, D. W. Balogh, L. A. Paquette, J. Am. Chem. pressure chemical ionization mass spectrometry measurements; length similar to the calculated HOMO-LUMO Soc. 104, 4503–4504 (1982). H. Takiguchi (The University of Tokyo) for microscopic Raman measurements; T. Tajima (Shimadzu) for solid-state UV-vis transition (278 nm), which can be caused by 3. H. W. Kroto, J. R. Heath, S. C. O’Brien, R. F. Curl, R. E. Smalley, absorption measurements and analysis; K. Oshida (The University Nature 318, 162–163 (1985). of Tokyo) for melting point measurements; and T. Aida (The vibronic coupling (fig. S11 and the supple- University of Tokyo), Y. Ito (The University of Tokyo), H. Sato 4. A. A. Popov, S. Yang, L. Dunsch, Chem. Rev. 113, 5989–6113 (RIKEN), K. Aikawa (The University of Tokyo), Y. Ishibashi mentary materials). (2013). (AGC Inc.), S. Kusumoto (The University of Tokyo), Y. Takemasa The reduction potentials of 1 and 5 (1: −2.1 V; (The University of Tokyo), Y. Nishibayashi (The University of 5. R. J. Cross, M. Saunders, H. Prinzbach, Org. Lett. 1, 1479–1481 Tokyo), S. Kuriyama (The University of Tokyo), H. Shinokubo 5: −2.5 V) were experimentally determined (1999). (Nagoya University), N. Fukui (Nagoya University), and using differential pulse voltammetry (DPV) Y. Segawa (Institute for Molecular Science) for assistance with 6. K. Komatsu, M. Murata, Y. Murata, Science 307, 238–240 experiments and fruitful discussions. We also thank U. Mayer at measurements (fig. S14). The LUMO levels of (2005). www.mayerscientificediting.com for proofreading. M.S. is grateful 1 and 5 were estimated to be −2.7 and −2.4 eV, to the Program for Leading Graduate Schools (MERIT-WINGS) respectively, corrected by comparison with 7. K. Kurotobi, Y. Murata, Science 333, 613–616 (2011). and to the JSPS for a fellowship for young scientists. A part of the HOMO level of ferrocene (−4.8 eV relative 8. K. K. Irikura, J. Phys. Chem. A 112, 983–988 (2008). this work was conducted in the Research Hub for Advanced to the vacuum level). Both values are in good 9. O. V. Boltalina, in New Fluorinated Carbons: Fundamentals and Nano Characterization (The University of Tokyo), which is supported by the Ministry of Education, Culture, Sports, Science agreement with the corresponding calculated Applications, O. V. Boltalina, T. Nakajima, Eds. (Elsevier, 2017), and Technology (MEXT; Japan), and in the One-stop Sharing values (1: −2.8 V eV; 5: −2.2 eV). In cyclic chap. 1. Facility Center for Future Drug Discoveries in the Graduate School voltammetry (CV) measurements, both 1 and 10. Y.-F. Wang et al., Theor. Chem. Acc. 127, 641–650 (2010). of Pharmaceutical Sciences (The University of Tokyo). Part of 5 showed an irreversible reduction wave in 11. F. Wahl et al., Chem. Eur. J. 12, 6255–6267 (2006). the theoretical calculations was carried out at the Research the range of −2.1 to −2.9 V (fig. S15), which 12. J. H. Holloway et al., J. Chem. Soc. Chem. Commun. 1991, Center for Computational Science, Okazaki, Japan. Funding: likely indicates that 1 and 5 accept an elec- 966–969 (1991). This work was financially supported by JSPS KAKENHI tron and that the resulting radical anions are 13. R. Taylor, J. Fluor. Chem. 125, 359–368 (2004). grants JP19K15532 (M.A.), JP21K14608 (M.A.), JP21J21713 14. J. Jia, C. Liu, H.-S. Wu, P. von Ragué Schleyer, H. Jiao, J. Phys. (M.S.), JP20H05839 (M.H.), and JP21H05385 (M.A.); the Toyota unstable under the applied conditions and Chem. C 113, 8077–8084 (2009). Riken Scholar Program; the Hattori Hokokai Foundation; the 15. R. Taylor et al., Nature 355, 27–28 (1992). Iketani Science and Technology Foundation; and the Koyanagi decompose. 16. T. Kato, T. Yamabe, J. Chem. Phys. 120, 1006–1016 Foundation (M.A.). AGC Inc. is acknowledged for supporting (2004). this work. Author contributions: M.A. conceived the concept To confirm the structure of the radical anion 17. E. W. Della, N. J. Head, J. Org. Chem. 60, 5303–5313 and directed the project. M.S., M.A., and Y.Y. conducted (1995). experiments. M.S. synthesized 1. M.A. and M.S. performed the generated upon accepting an electron into the 18. H. Irngartinger, S. Strack, J. Am. Chem. Soc. 120, 5818–5819 crystallographic studies and the DFT calculations. K.K. LUMO of 1, derived from the eight C–F s* or- (1998). performed the matrix-isolation ESR measurements. M.H. bitals, we used matrix-isolation electron spin 19. H. H. Freedman, D. R. Petersen, J. Am. Chem. Soc. 84, performed Frank-Condon and Herzberg-Teller analysis. All 2837–2838 (1962). authors analyzed and discussed the results. M.S. and M.A. resonance (ESR) spectroscopy in combination 20. K. F. Biegasiewicz, J. R. Griffiths, G. P. Savage, J. Tsanaktsidis, prepared the manuscript with feedback from all authors. K.N. with g-ray radiolysis at low temperatures. This R. Priefer, Chem. Rev. 115, 6719–6745 (2015). and T.O. supervised the project. Competing interests: The technique has already been used to observe un- 21. M. J. Gerace, D. M. Lemal, H. Ertl, J. Am. Chem. Soc. 97, authors declare that they have no competing interests. Data 5584–5586 (1975). and materials availability: The x-ray data have been deposited stable radical anions of perfluorocycloalkanes 22. A. C. BarefootIII, W. D. Saunders, J. M. Buzby, M. W. Grayston, in the Cambridge Crystallographic Data Centre (CCDC) under (36). As shown in Fig. 4C, the ESR spectrum D. M. Lemal, J. Org. Chem. 45, 4292–4295 (1980). reference numbers 2144167 to 2144171. All other data are of 1¥− in a hexamethylethane (HME) matrix at 23. J. L. Adcock, H. Zhang, J. Org. Chem. 61, 1975–1977 presented in the main text or the supplementary materials. 77 K is well resolved. Although the innermost (1996). License information: Copyright © 2022 the authors, some 24. T. Okazoe, J. Fluor. Chem. 174, 120–131 (2015). rights reserved; exclusive licensee American Association for the peaks overlap with the intense signals of the 25. G. Sandford, J. Fluor. Chem. 128, 90–104 (2007). Advancement of Science. No claim to original US government 26. K. Murata, S.-Z. Wang, Y. Morizawa, K. Oharu, J. Fluor. Chem. works. https://www.science.org/about/science-licenses- HME matrix radicals in the congested central 127, 1125–1129 (2006). journal-article-reuse 27. S. Berski, A. J. Gordon, Z. Latajka, J. Phys. Chem. A 118, region, the other eight peak bands were ob- 4147–4156 (2014). SUPPLEMENTARY MATERIALS served in accordance with the total nuclear science.org/doi/10.1126/science.abq0516 spin quantum number, MI = ±4, ±3, ±2, and Materials and Methods ±1. Because of the large a(19F) hyperfine Supplementary Text Figs. S1 to S67 Tables S1 to S6 References (39–56) Submitted 15 March 2022; accepted 7 June 2022 10.1126/science.abq0516 SCIENCE science.org 12 AUGUST 2022 • VOL 377 ISSUE 6607 759

RESEARCH | REPORTS SPEECH EVOLUTION (CT) to examine larynges from 25 genera and 43 species of primates (table S1). All non- Evolutionary loss of complexity in human vocal human species possess a vocal membrane anatomy as an adaptation for speech and often exhibit a shallow sulcus separating the vocal membrane from the vocal fold (Fig. Takeshi Nishimura1,2*, Isao T. Tokuda3, Shigehiro Miyachi1,2, Jacob C. Dunn4,5,6, Christian T. Herbst1,6, 1C). Platyrrhines usually have a tall vocal Kazuyoshi Ishimura3, Akihisa Kaneko1,2, Yuki Kinoshita1,2, Hiroki Koda1 , Jaap P. P. Saers5, membrane (Fig. 1, D to G, and fig. S1, A to Hirohiko Imai7, Tetsuya Matsuda7, Ole Næsbye Larsen8, Uwe Jürgens9, Hideki Hirabayashi10, G), whereas cercopithecids have a shorter lip- Shozo Kojima1, W. Tecumseh Fitch6,11* like membrane (Fig. 1, D and H to J, and fig. S1, H to L). This feature can vary even within a Human speech production obeys the same acoustic principles as vocal production in other animals single species in hominoids (Fig. 1, C, D, and but has distinctive features: A stable vocal source is filtered by rapidly changing formant frequencies. K to M; and fig. S1, M to P). In the gibbons To understand speech evolution, we examined a wide range of primates, combining observations of and siamang (hylobatids), the vocal membrane phonation with mathematical modeling. We found that source stability relies upon simplifications in extends from the lateral wall of the laryngeal laryngeal anatomy, specifically the loss of air sacs and vocal membranes. We conclude that the cavity and has become disconnected from the evolutionary loss of vocal membranes allows human speech to mostly avoid the spontaneous nonlinear vocal fold (Fig. 1K and fig. S1, M to O). phenomena and acoustic chaos common in other primate vocalizations. This loss allows our larynx to produce stable, harmonic-rich phonation, ideally highlighting formant changes that convey most Given these anatomical data, the most par- phonetic information. Paradoxically, the increased complexity of human spoken language thus followed simonious evolutionary conclusion is that the simplification of our laryngeal anatomy. vocal membrane is an ancestral primate fea- ture, and that this feature was lost in the hu- S peech is the dominant mode of human nal provides the backdrop upon which rapid man lineage to yield the rounded vocal fold linguistic expression, and for most of modulations of our vocal tract filter yield a typical of humans (Fig. 1B). Thus, we argue that dynamic pattern of formant frequencies con- the absence of a vocal membrane in the human our evolutionary history, until the emer- veying phonetic information. Together, under larynx is an evolutionarily derived feature. fine neural control, this human source–filter gence of signed and written languages, system provides the high-bandwidth commu- Turning to function, the anatomical features nicative signal needed to rapidly encode com- of the vocal membrane suggest that it plays a speech provided the sole communicative plex linguistic information. role in phonation. When the glottis (the air space between the vocal folds) is closed for modality for language. Speech production is Speech-related specializations of the human phonation by adducting the arytenoid carti- vocal tract, including the descent of the tongue lages, both the vocal folds and membranes based on the same acoustic and physiological root into the pharynx, are well documented move toward the glottal midline, implying (10–12), and their acoustic effects well under- that the vocal membranes should also vibrate principles as vocal production in other terres- stood (5, 8–11), but evolutionary changes in during phonation in primates. trial vertebrates (1–4) but nonetheless pos- our larynx have been relatively neglected. First, sesses distinctive attributes. First, our vocal humans lost the laryngeal air sacs seen in Supporting this hypothesis, vocal membrane other great apes (6, 13), and which were prob- vibration was documented in vivo during pho- source (produced by laryngeal vocal fold oscil- ably still present in Australopithecus (14). Sec- nation in a chimpanzee, the closest phylo- ond, we show here that vocal membranes (also genetic relative of humans. Hayama et al. lations modulating air flow from the lungs; known as “vocal lips”)—thin upward projec- examined reflex glottal closure using a trans- tions of the vocal folds (Fig. 1C)—are typical nasal fiberscope in an adult male chimpanzee Fig. 1, A and B) is uncharacteristically stable, laryngeal features in primates but were lost under anesthesia (15). Serendipitously, these in humans (Fig. 1B). Using a combination of video data also documented phonation during and nonpathological adult speech completely multidisciplinary methods, we show that vocal 375 grunts or growls as the chimpanzee was membranes increase nonlinearities, yielding awakening from anesthesia. We reexamined avoids the nonlinear phenomena and bifurca- vocal instability. This leads to the surprising these data and found that the vocal mem- tions to chaos—i.e., irregular oscillations and conclusion that the increased stability of hu- branes always vibrate during vocalization, abrupt frequency transitions—commonly seen man phonation results from an evolutionary colliding to close the glottis (Fig. 2A and in most other mammals (5, 6). This stability, loss of anatomical complexity. Although fossil movie S1). These observations suggest that, indicators of vocal fold anatomy are unavail- in chimpanzees, the vocal membrane plays a combined with enhanced neural control of able, our comparative data indicate that this central role in phonation even during low- laryngeal muscles (7), yields a highly reliable simplification through loss must have occurred frequency calls and is not simply an accessory fundamental frequency (fo or “pitch”) and rich since our divergence from chimpanzees roughly feature of the vocal folds subserving high-pitch array of harmonics (5, 8, 9). The predictable, 6 million years ago. vocalization as previously thought (5, 16, 17). broadband acoustic energy in this source sig- We show that a vocal membrane is a key Because such in vivo evidence is very chal- anatomical feature shared by all nonhuman lenging to obtain in chimpanzees, we next 1Primate Research Institute, Kyoto University, Inuyama, Aichi anthropoid primates (“primates” hereafter), in- quantitatively investigated the role of vocal 484-8506, Japan. 2Center for the Evolutionary Origins of Human cluding hominoids or apes, cercopithecids membranes in an ex vivo setting, using ex- Behavior, Kyoto University, Inuyama, Aichi 484-8506, Japan. or Old World monkeys, and platyrrhines or cised larynges from three chimpanzees. Each 3Department of Mechanical Engineering, Ritsumeikan University, New World monkeys (see Fig. 1D for phyletic larynx was mounted on a vertical tube supplying Kusatsu, Shiga 525-8577, Japan. 4Behavioural Ecology Research relationships). We used magnetic resonance airflow, and the adductory gestures observed Group, School of Life Science, Anglia Ruskin University, imaging (MRI) and computed tomography in vivo were reproduced by positioning the Cambridge CB1 1PT, UK. 5Department of Archaeology, arytenoid cartilages with adjustable prongs. University of Cambridge, Cambridge CB2 3DZ, UK. 6Department Adductory conditions were further modified of Behavioral and Cognitive Biology, University of Vienna, 1030 by pulling the thyroid cartilage antero-inferiorly. Vienna, Austria. 7Department of Systems Science, Graduate Vocal membrane and fold vibrations driven by School of Informatics, Kyoto University, Sakyo, Kyoto 606-8501, Japan. 8Department of Biology, University of Southern Denmark, DK-5230 Odense M, Denmark. 9Section of Neurobiology, German Primate Center, D-37077 Göttingen, Germany 10Dokkyo Medical University, Mibu, Tochigi 321-0293, Japan. 11Cognitive Science Hub, University of Vienna, Vienna, Austria. *Corresponding author. Email: [email protected] (T.N.); [email protected] (W.T.F.) Present address: Department of Life Science, Graduate School of Arts and Sciences, The University of Tokyo, Meguro, Tokyo 153-8902, Japan. 760 12 AUGUST 2022 • VOL 377 ISSUE 6607 science.org SCIENCE

RESEARCH | REPORTS varying air pressure were documented by using Fig. 1. Vocal anatomy and its phylogeny in anthropoids. (A) The cycle of vocal fold vibration in wave high-speed video, and phonatory dynamics fashion in humans. (B) Frontal section (gradation inversion, left) and a corresponding line drawing (middle) were simultaneously documented by using at the level of the dashed line on an excised and formalin-fixed larynx (right) in humans (courtesy of K. Sato). time-synchronized acoustic and electroglot- (C) Frontal MRI scan (left) and a corresponding line drawing (middle) at the level of the dashed line on vocal tographic recordings. anatomy reconstructed from MRI scans (right) for a chimpanzee, Pan troglodytes. Blue indicates the underlying arytenoid cartilage and green the thyroarytenoid muscle. (D) Variation and phylogeny of the vocal The vocal membranes always participated in membrane and (E to M) MRIÐCT frontal scans in the species labeled on (D). Transverse scale bar, 5 mm; vibration in these chimpanzee ex vivo experi- vertical scale bar, 2.5 mm. Key: ab, body of the arytenoid cartilage; s, sulcus; ta, thyroarytenoid muscle; vef, ments (Fig. 2, B to G, and movie S2, A to C) and ventricular fold; vf, vocal fold; vm, vocal membrane; and vp, vocal process of the arytenoid cartilage. The oscillated in a wave-like fashion in the absence asterisk (*) indicates a posterior commissure of the vocal membrane and fold. of vocal fold vibration (Fig. 2, B and C, and movie S2A). The vocal membranes also vi- quired further investigation of vocal membrane and S8), and also ex vivo in our chimpanzee brated separately from—but in phase with—the function. Although many primate species can (Fig. 2, F and G, and movie S2C) and rhesus vocal folds (Fig. 2, D and E, and movie S2B). In produce clear harmonically structured calls macaque data (fig. S2, I and J, and movie such cases, the vocal membranes usually col- [e.g., marmoset phees (17), macaque coos (18), S6B). The chaotic episodes were associated lided, but the vocal folds rarely did so. Only or chimpanzee hoos (19)], even these call types with highly irregular vibrations and collisions when very firmly adducting and lowering often bifurcate to subharmonics or chaos at of the vocal membranes, superimposed upon the arytenoid vocal process did we observe higher intensities (6, 17). Our findings below nearly periodic vibrations of the vocal folds vocal fold collision resembling that of humans. indicate that the primate vocal membrane in chimpanzees ex vivo (Fig. 3, B and C, and In such cases, there was sometimes a phase- plays an important acoustic role by increasing movies S2D and S3E) and macaques in vivo delayed oscillatory pattern, in which the vocal susceptibility to nonlinear phenomena. (fig. S3, A to C, and movie S5, C and D). This folds first collided and the vocal membranes suggests that vocal membrane vibrations play followed (Fig. 2, F and G, and movie S2C). Subharmonics and chaos were observed a crucial role in generating the high propor- Thus, adducting the arytenoid always resulted in vivo in both our rhesus macaque (fig. S3, A tion of subharmonics and chaos empirically in vocal membrane vibration and often colli- to C, and movie S5, C and D) and squirrel observed in nonhuman primate vocalization. sion. By contrast, vibrations of the vocal folds monkey data (fig. S3, D and E, and movies S7 alone, without the vocal membrane, were never observed in these experiments. In sharp con- trast, vocal fold vibrations always play a central role in phonation in human speech. Further in vivo experiments with two rhesus macaques (cercopithecids) and two squirrel monkeys (platyrrhines) corroborated the ob- servations in chimpanzees. We successfully induced vocalization by electrical stimulation to the periaqueductal gray and surrounding areas of the midbrain in these monkeys under anesthesia. In both species, high-speed video recordings showed that the vocal membranes always vibrated and, in many cases, collided during phonation (Fig. 2, H and I, and movies S3 to S4). In macaques, as in chimpanzees, the vocal folds did not always vibrate (fig. S2, A and B, and movie S5A), but the second subject did show simultaneous vibrations of both vocal membranes and folds, either in phase (fig. S2, C and D, and movie S5B) or phase delayed (fig. S2, E and F, and movie S5C). Further ex vivo experiments with six macaque larynges repro- duced the same patterns (fig. S2, G to L, and movie S6). Humanlike vibration of the vocal folds alone was not observed in macaques ex vivo. These physiological data, combined with previous work, provide clear empirical evidence that vocal membranes constitute the predominant vocal source generator in pri- mates, always vibrating and typically collid- ing, whereas the contributions of the vocal folds are limited or even absent. Taken together, the evolutionary loss of the vocal membrane transformed the predomi- nant source generator from the vocal mem- branes to the vocal folds in the human lineage. However, understanding the adaptive impor- tance of the derived human condition re- SCIENCE science.org 12 AUGUST 2022 • VOL 377 ISSUE 6607 761

RESEARCH | REPORTS Fig. 2. Vocal membrane and vocal fold vibrations. (A) Images extracted from video recordings in vivo in a chimpanzee. From left to right: The vocal membranes and folds are maximally separated, the glottis is closing, and the vocal membranes close the glottal space and vibrate. (B to G) Vibration of the vocal membrane and vocal fold ex vivo in chimpanzees. (B, D, and F) A top-view image of the glottis extracted from video recordings with the dotted line for creating the kymograph (left) and kymograph showing a time series of glottal vibration (right). (C, E, and G) Electroglottograph (upper panel) and acoustic (lower panel) signals. (B and C) Only the vocal membranes vibrate, whereas the vocal folds do not; and (D and E) vocal membrane and fold both vibrate in-phase and (F and G) out-of-phase. (H and I) Images extracted from video recordings in vivo in a rhesus macaque, Macaca mulatta (H) and a common squirrel monkey, Saimiri sciureus (I). See Fig. 1 for the anatomical key. We used a mathematical model to deter- empirical data from marmosets (17), suggest- freedom to laryngeal biomechanics. This dy- mine the acoustic effects of vocal membrane ing that the vocal membrane provides im- namically generates multiple routes to voice vibration in primates compared with vocal proved efficiency in phonation and/or allows instabilities, as seen in many coupled non- fold vibration alone, as typifies humans. In our louder and higher-frequency vocalizations in linear systems (22–24). Thus, the simulations model, the vocal fold is represented by two nonhuman primates. confirmed our empirical observations, indicat- masses coupled with springs (3, 16, 20–22) ing that the primate vocal membrane, com- (Fig. 3D). The vocal membrane is modeled as In bifurcation diagrams, we showed that the bined with vocal fold vibration, destabilizes an additional reed-like plate that can vibrate nonhuman model gives rise to various non- the vocal source. Their interaction can spon- independently, attached to the upper mass linear phenomena, including subharmonics taneously lead to subharmonics and chaos in (21) (Fig. 3D). We first confirmed that all and chaos, with increasing adduction (Fig. 3E), response to simple linear variations in glottal three empirically observed vibration patterns reproducing the chaotic episode observed in and respiratory parameters. The loss of the can be reproduced by the numerical simula- chimpanzees ex vivo (Fig. 3, B and C, and vocal membrane in humans therefore reduces tions in the nonhuman model with vocal movies S2D and S3E). Similar bifurcations the risk of contaminating the stable vocal fold membranes (fig. S3F). Next, we generated bi- were observed with increasing subglottal pres- oscillations used in human speech or singing furcation diagrams by adducting the vocal sure in the nonhuman model (fig. S3G). This with chaotic irregularities and noise. membranes and folds in the model (Fig. 3, E contrasts with the situation in our human and F). Notably, the phonation threshold pres- model, which generated only stable periodic In summary, the primate vocal membrane sure observed at 0.14 kPa in the nonhuman vibrations regardless of increasing adduction interacts with the vocal fold beneath it to ef- model was much lower than that of the hu- (Fig. 3F) and pressure (fig. S3H). ficiently generate phonation but simultane- man model lacking vocal membranes (0.19 kPa). ously readily generates nonlinear phenomena This is consistent with previous theoretical The increased instability seen in the non- that destabilize the vocal source. The resulting models of a fixed vocal membrane (16) and human model is not surprising, because the spontaneous phenomena are unlikely to be addition of mechanically coupled vocal mem- under the animal’s volitional control, but this branes contributes additional degrees of 762 12 AUGUST 2022 • VOL 377 ISSUE 6607 science.org SCIENCE

RESEARCH | REPORTS Fig. 3. Bifurcation observed and numerical simulation of mathematical model. (A) Sound spectrograph 19. P. Marler, R. Tenaza, in How Animals Communicate, T. A. Sebeok, showing bifurcation to chaos ex vivo in a chimpanzee. (B) Kymographs (left) and electroglottograph Ed. (Indiana Univ. Press, Bloomington, 1977), pp. 965–1033. signals (right) for a stable phonation at 1.2 s of (A); and (C) for the chaotic phonation at 2.7 s of (A). (D) Schematic illustration of the nonhuman model with a vocal membrane. (E) Bifurcation diagram created 20. K. Ishizaka, J. L. Flanagan, Bell Syst. Tech. J. 51, 1233–1268 by slowly adducting (increasing the contact area) for the nonhuman and (F) the human models. See Fig. 1 (1972). for the anatomical key. 21. J. Neubauer, “Nonlinear dynamics of the voice: Bifurcations more complex laryngeal anatomy enables non- neural control mechanisms, represented a cru- and mode analysis of complex spatio-temporal signals,” thesis, human primates to efficiently generate a di- cial evolutionary transition supporting the Humboldt University of Berlin (2004). verse set of complex vocal signals without the sophistication and complexity of human spo- need for correspondingly complex motor control. ken language. 22. I. T. Tokuda, Anthropol. Sci. 126, 35–41 (2018). 23. P. Berg, Y. Pomeau, C. Vidal, Order within Chaos: Towards a By contrast, despite its disadvantages for REFERENCES AND NOTES efficiency, we suggest that the stable low- Deterministic Approach to Turbulence (Wiley, 1986). pitched phonation resulting from the loss of 1. A. M. Taylor, D. Reby, J. Zool. 280, 221–236 (2010). 24. L. Glass, M. C. Mackey, From Clocks to Chaos: The Rhythms of the vocal membranes in human evolution en- 2. C. T. Herbst et al., Science 337, 595–599 (2012). hances the detectability and salience of for- 3. H. Koda, I. T. Tokuda, M. Wakita, T. Ito, T. Nishimura, Life (Princeton Univ. Press, 1988). mant frequencies, which carry most of the 25. M. J. Owren, D. Rendall, J. Comp. Psychol. 117, 380–390 (2003). phonetic information in human speech. Chaos J. Acoust. Soc. Am. 137, 3068–3076 (2015). 26. W. T. Fitch, J. B. Fritz, J. Acoust. Soc. Am. 120, 2132–2141 and subharmonics generate an irregular and 4. C. P. Elemans et al., Nat. Commun. 6, 8978 (2015). changing spectrum thought to obscure the 5. W. T. Fitch, Trends Cogn. Sci. 4, 258–267 (2000). (2006). conformation of the vocal tract filter (18, 25) 6. W. T. Fitch, J. Neubauer, H. Herzel, Anim. Behav. 63, 407–418 27. W. Mende, H. Herzel, K. Wermke, Phys. Lett. A 145, 418–424 by introducing spurious spectral peaks or “pseudo-formants” (26). Although nonlinear (2002). (1990). phenomena can occur in humans—e.g., in 7. U. Jürgens, Neurosci. Biobehav. Rev. 26, 235–258 28. H. Herzel, D. Berry, I. R. Titze, M. Saleh, J. Speech Hear. Res. screams and baby cries (27)—they are rarely observed during speech except in patholog- (2002). 37, 1008–1019 (1994). ical cases (28). We conclude that, in addition 8. I. R. Titze, Principles of Voice Production (Prentice Hall, 1994). to changes in vocal tract anatomy and neural 9. K. N. Stevens, Acoustic Phonetics (MIT Press, 1998). ACKNOWLEDGMENTS control, this previously unnoticed simplifi- 10. P. H. Lieberman, D. H. Klatt, W. H. Wilson, Science 164, cation of human laryngeal anatomy allowed We are grateful to the late Sugio Hayama for conceiving and our species to easily generate and regulate 1185–1187 (1969). conducting in vivo chimpanzee experiments and for establishing complex, high-bandwidth vocalizations under 11. P. Lieberman, Curr. Anthropol. 48, 39–66 (2007). the collection of excised larynges at the Japan Monkey Centre cognitive control in speech. Thus, seemingly 12. T. Nishimura, A. Mikami, J. Suzuki, T. Matsuzawa, Proc. Natl. (JMC). We thank the JMC, Nagoya Higashiyama Zoo, Fukuoka City paradoxically, a loss of complexity in human Zoological Garden, Kagoshima City Hirakawa Zoological Park, laryngeal anatomy, coupled with enhanced Acad. Sci. U.S.A. 100, 6930–6933 (2003). Yokohama Zoological Gardens, and the Wildlife Research Center, 13. S. Hayama, J. Anthrol. Soc. Jpn. 78, 274–298 (1970). the Laboratory of Physical Anthropology, and the Primate Research 14. Z. Alemseged et al., Nature 443, 296–301 (2006). Institute (PRI) of Kyoto University for providing samples. We thank 15. S. Hayama, in Primatology Today, A. Ehara, T. Kimura, T. Udono, Y. Shintaku, K. Matsui, H. Imai, S. Tsubouchi, and A. Kataoka for preparing or scanning the specimens; Y. Nomura and Z. Chen for O. Takenaka, M. Iwamoto, Eds. (Elsevier, 1991), pp. 491–492. examining laryngeal anatomy; K. Migimatsu, T. Matsumoto, 16. P. Mergell, W. T. Fitch, H. Herzel, J. Acoust. Soc. Am. 105, M. Kanaya, and R. Miyazaki for conducting ex vivo experiments; C. R. Larson for providing suggestions for in vivo experiments of 2020–2028 (1999). macaques; and T. Miyabe-Nishiwaki, Y. Sawada, B. Bach Andersen, and 17. Y. S. Zhang, D. Y. Takahashi, D. A. Liao, A. A. Ghazanfar, present and past staff of the PRI and the German Primate Center for help conducting in-vivo experiments and/or daily care of subjects. C. P. H. Elemans, Nat. Commun. 10, 4592 (2019). We appreciate the photograph courtesy of K. Sato at the Kurume 18. D. Rendall, M. J. Owren, P. S. Rodman, J. Acoust. Soc. Am. 103, University. This research was in part supported by the Kyoto University- University of Vienna Strategic Partnership Program. Funding: Japan 602–614 (1998). Society for the Promotion of Science grants 16H04848 and 19H01002 (T.N.), 17H06313 and 20K11875 (I.T.T.), and 18H03503 (H.K.) Research Units for Exploring Future Horizons through the Kyoto University Research Coordination Alliance (T.N., C.T.H.) Rhinology and Laryngology Research Fund UK (J.C.D.) Ministry of Education, Culture, Sports, Science and Technology Grant-in-aid for Scientific Research on Innovative Areas no. 4903 (Evolinguistics) grant 17H06380 (H.K.) Ministry of Education, Culture, Sports, Science and Technology Grant-in-aid for Scientific Research on Innovative Areas no. 4903 (Evolinguistics) (W.T.F.) Austrian Science Fund DK Grant “Cognition & Communication 2” grant W1262-B29 (W.T.F.). Author contributions: Conceptualization: T.N., W.T.F. Methodology: T.N., C.T.H., I.T.T., S.M., J.C.D., U.J., W.T.F., O.N.L. Formal analysis: T.N., I.T.T., W.T.F. Investigation: T.N., W.T.F., I.T.T., S.M., J.C.D., C.T.H., K.I., A.K., Y.K., H.K., J.P.P.S., H.I., T.M., O.N.L., U.J., H.H., S.K. Resources: T.N., A.K., U.J., S.K., O.N.L. Writing - original draft: T.N., I.T.T., S.M., J.C.D., C.T.H., W.T.F. Writing - review and editing: T.N., I.T.T., W.T.F. Project administration: T.N. Funding acquisition: T.N., I.T.T., J.C.D., H.K., W.T.F. Competing interests: The authors declare that they have no competing interests. Data and materials availability: Laryngeal specimens from the subjects named Baran, Keiko, Wilie, and Yuri were available from the zoos under a material transfer agreement with the Kyoto University. Some specimens were provided through the Collaborative Research program of the JMC (no. 2018017) and the Great Ape Information Network under the National BioResource Project (NBRP) of Japan. All data are available in the main text or the supplementary materials. License information: Copyright © 2022 the authors, some rights reserved; exclusive licensee American Association for the Advancement of Science. No claim to original US government works. https://www. science.org/about/science-licenses-journal-article-reuse SUPPLEMENTARY MATERIALS science.org/doi/10.1126/science.abm1574 Materials and Methods Supplementary Text Figs. S1 to S4 Table S1 References (29–35) Movies S1 to S8 Submitted 30 August 2021; accepted 30 June 2022 10.1126/science.abm1574 SCIENCE science.org 12 AUGUST 2022 • VOL 377 ISSUE 6607 763

RESEARCH | REPORTS INSECT MIGRATION n = 8) and were then followed for a mini- mum of 1 hour and up to 3.65 hours (2.5 ± Individual tracking reveals long-distance flight-path 0.30 hours, n = 8; table S1). The moths were control in a nocturnally migrating moth followed for a mean distance of 62.7 ± 6.7 km (n = 8) and up to 89.6 km (Fig. 1C and tables Myles H. M. Menz1,2,3*, Martina Scacco1,3, Hans-Martin Bürki-Spycher4, Hannah J. Williams1,3, S1 and S2), the longest distance over which Don R. Reynolds5,6, Jason W. Chapman7,8,9*, Martin Wikelski1,3,10* any insect has been continuously tracked in the field. The overall migration direction was Each year, trillions of insects make long-range seasonal migrations. These movements are relatively well toward the south-southwest [Rayleigh test: understood at a population level, but how individual insects achieve them remains elusive. Behavioral 208.70° ± 0.42° (mean ± SD), mean vector responses to conditions en route are little studied, primarily owing to the challenges of tracking length (r) = 0.917, P ≤ 0.001, n = 14; Fig. 1D]. individual insects. Using a light aircraft and individual radio tracking, we show that nocturnally migrating This track direction is very similar to the deathÕs-head hawkmoths maintain control of their flight trajectories over long distances. The moths preferred headings of a range of migratory did not just fly with favorable tailwinds; during a given night, they also adjusted for head and crosswinds insects (moths, butterflies, and hoverflies) that to precisely hold course. This behavior indicates that the moths use a sophisticated internal compass were observed with radar in Western Europe to maintain seasonally beneficial migratory trajectories independent of wind conditions, illuminating how (2, 10, 19, 20), including hawkmoths (10), all of insects traverse long distances to take advantage of seasonal resources. which likely follow a similar western route to the Mediterranean or northwest Africa. I nsect migration takes place on an enor- and modify them with respect to ambient mous scale, with trillions of individuals wind conditions. However, owing to the meth- We obtained detailed tracks for seven of odological constraints of tracking such small these moths, each with three or more locations performing bidirectional seasonal move- animals over long distances at night (15), in a single night (table S2). Moths traveled individual moths have never been tracked with a mean ground speed of 9.4 ± 0.4 m/s ments that have important impacts on throughout their migration, and so the capa- (33.8 km/hour; n = 99 segments; Fig. 1E) and a bility of these migrants to maintain straight maximum recorded ground speed of 19.4 m/s ecosystem function and provision of es- flight paths, over long distances and in sea- (69.7 km/hour). The mean ground speed re- sential services (1–5). However, the naviga- sonally beneficial directions, is unknown. corded (Fig. 1E) is consistent with what we tional mechanisms and behavioral strategies expect the upper limit of self-powered flight in We used animal-borne radio telemetry to A. atropos to be (21), suggesting that moths used by night-flying migrants, especially larger record complete tracks of individually tagged modulated their self-powered airspeed and/or moths over a full night during autumn mi- received relatively modest wind assistance. Al- nocturnal lepidopterans (macromoths), dur- gration, within the context of the fine-scale though there was variation in individual mi- wind fields experienced as they migrated south- gration direction, all moths maintained straight ing these long-range journeys have been un- ward through the Alps of central Europe. tracks (straightness index: mean = 0.95, range = Our study species, the death’s-head hawkmoth 0.80 to 0.99, n = 7; Fig. 1C and table S2) along known for more than 100 years. (Acherontia atropos, Sphingidae; Fig. 1A), is their entire flight paths, which lasted many Europe’s largest lepidopteran, with a rich folk- tens of kilometers, despite being subjected The view in the first half of the 20th century, lore that stems from its sinister skull-like to winds of varying strength and direction thoracic markings, unusual habit of raiding throughout their course (Fig. 2). Two of the promoted by C. B. Williams, was that migrant beehives to steal honey, and startling acoustic seven moths evidently crossed the Alps during capabilities (16, 17). A. atropos is a long- a single night, because they were relocated moths controlled their movement direction distance Afro-Palearctic migrant, arriving to south of the Alps during searches early the breed in Europe north of the Alps each spring. following morning. Their locations were con- irrespective of the wind and maintained The subsequent generation returns south the sistent with their individual trajectories re- straight flightpaths over long distances (6, 7). following autumn to winter-breeding regions corded the preceding night, suggesting that in the Mediterranean Basin and likely also they had maintained straight tracks even Empirical evidence of persistent, self-directed sub-Saharan Africa (16, 17), covering a distance while transiting the Alps [covering distances of up to 4000 km. The moths are extremely of 173.9 and 161.8 km from the release point tracks was lacking, however, and by the sec- large for flying insects, weighing up to 3.5 g (Fig. 1C and table S1)]. [2.65 ± 0.15 g (mean ± SE), n = 14] and capable ond half of the 20th century, C. G. Johnson of carrying tiny very high frequency (VHF) To answer the question of how moths are radio transmitters. We used a light aircraft able to maintain straight tracks relative to the and L. R. Taylor downplayed the importance (Cessna 172) to track hawkmoths fitted with ground while exposed to varying winds, we transmitters (Fig. 1A) and recorded precise calculated the distribution of the angle of de- of orientation behavior and emphasized the (±150 m) GPS locations from the aircraft (18) viation, b (the difference between the track at regular intervals throughout their migra- and the downwind direction), to determine role of wind in determining migratory trajec- tion (every 5 to 15 min, when possible). the extent to which the self-powered heading tories (8, 9). The modern view has swung back influenced the trajectory (22). The analysis We recorded nocturnal migratory flights of revealed that moths used three distinct be- again, because radar observations of free- 14 moths, with eight at high spatiotemporal havioral strategies, which resulted in the flight flying migrants (10–12) and experimental ma- frequency, as they migrated toward the Med- paths of the moths grouping into three direc- nipulation of tethered individuals (13, 14) have iterranean (Fig. 1, B and C, and table S1). tional clusters (Fig. 1, C and D). These clusters Moths initiated migration at a similar time appeared to be partly determined by the am- both clearly demonstrated that nocturnally after sunset (62 ± 4.9 min, range 42 to 81 min, bient wind conditions experienced along the flight path (Fig. 2) and partly by the to- migrating moths can select adaptive headings pography of the landscape (Fig. 1C). 1Department of Migration, Max Planck Institute of Animal Behavior, 78315 Radolfzell, Germany. 2College of Science and Engineering, James Cook University, Townsville, QLD 4811, Australia. 3Department of Biology, University of Konstanz, 78464 Konstanz, Germany. 4Independent Researcher, Promenadenstrasse 2, 3076 Worb, Switzerland. 5Natural Resources Institute, University of Greenwich, Chatham, Kent ME4 4TB, UK. 6Rothamsted Research, Harpenden, Hertfordshire AL5 2JQ, UK. 7Centre for Ecology and Conservation, University of Exeter, Penryn, Cornwall TR10 9FE, UK. 8Environment and Sustainability Institute, University of Exeter, Penryn, Cornwall TR10 9FE, UK. 9Department of Entomology, College of Plant Protection, Nanjing Agricultural University, Nanjing 210095, China. 10Centre for the Advanced Study of Collective Behaviour, University of Konstanz, 78464 Konstanz, Germany. *Corresponding author. Email: [email protected] (M.H.M.M.); [email protected] (M.W.); j.chapman2@exeter. ac.uk (J.W.C.) 764 12 AUGUST 2022 • VOL 377 ISSUE 6607 science.org SCIENCE

RESEARCH | REPORTS Fig. 1. Individual tracking reveals consistently straight flight paths in migrating hawkmoths. (A) Death’s- The first strategy was used under opposing head hawkmoth (A. atropos) showing the characteristic skull-like marking (left) and fitted with a miniaturized wind directions and resulted in moths taking VHF radio transmitter (weighing 240 mg) (right). (B) Map of Europe with the tracking region outlined [magnified the most direct route to the wintering grounds view shown in (C)]; the release site was ~50 km north of the Swiss Alps. (C) Nighttime tracks of migrating by maintaining a constant southward track hawkmoths showing persistently straight trajectories throughout a night’s flight. Solid lines indicate moths that (Fig. 2, A and B). Under this strategy, the were tracked continuously throughout a night, and dashed lines indicate presumed tracks of moths to their moths continuously adjusted their headings relocation position during searches in the following day(s). Colors represent different individuals and are consistent so that distributions of b had 95% confidence between figures. (D) Mean track directions of the 14 moths that demonstrated migratory behavior. Each point intervals (CIs) that overlapped 180° and had a represents the track direction of an individual moth. The arrow indicates the overall mean direction (208.70°), and mean b close to that value (Fig. 3A and table S2), arrow length indicates the directedness (r = 0.917). (E) Frequency distribution of track segments in each of resulting in more-or-less upwind flight (Fig. 4, A 10 ground-speed categories (each bar represents a range of 2 m/s) for the seven moths that were to D). Examination of ground speeds and wind tracked continuously on migration (n = 99 segments). The dashed line indicates the mean of 9.4 m/s. speeds along the track (figs. S1 and S2) (21) [Photo credit: Christian Ziegler] indicated that moths that used this orientation strategy must have flown close to the ground (50 m or lower), that is, within their “flight boundary layer” [the lower-most layer of the atmosphere within which the insect’s self- powered flight speed exceeds the wind speed, allowing control of their trajectory (4, 22)]. The second and third strategies were both used under favorable wind directions (i.e., oc- casions when southward flight would expose moths to some degree of tailwind assistance). We predicted that moths using tailwind assist- ance would fly in the layer where winds were fastest, as previously observed in studies of noctuid moths (10, 19). However, examination of ground speeds and airspeeds on these oc- casions indicated that hawkmoths that used these strategies flew about 300 m above the ground, considerably lower than the wind- speed maxima available (figs. S1 and S2) but high enough to receive some wind assis- tance (Fig. 4, E to G). Under these conditions, moths appear to balance speed with direction, as seen in other migrant moths (10). The sec- ond orientation strategy involved flying rela- tively close to the south-westward downwind direction (Fig. 2, C and D), but individuals modified their heading to achieve a straight trajectory lying somewhat further south of the strongest wind (as supported by values of b around −30° to −50° and for which the 95% CIs do not overlap with 0°; Fig. 3B and table S2). The final orientation strategy, which was used by a single individual (moth 5), involved flying directly downwind (as indicated by the 95% CI of b overlapping 0°; Fig. 3C and table S2), resulting in a track toward the west-southwest (Fig. 2D) with a higher ground speed than any other moth (Fig. 4H). In general, there was a negative relation- ship between airspeed and wind assistance, with airspeed increasing in headwinds and decreasing in tailwinds (Fig. 4). Furthermore, median ground speed was relatively similar across the orientation strategies (Fig. 4H). Thus, moths modulated their ground speed by varying their self-powered flight vector under different wind conditions to achieve a pre- ferred ground speed, similar to that docu- mented in many insects (23), which may be beneficial in the trade-off between energy consumption and travel speed (22). SCIENCE science.org 12 AUGUST 2022 • VOL 377 ISSUE 6607 765

RESEARCH | REPORTS Fig. 2. Migrating hawkmoths continuously compensate for wind to maintain southwest and thus skirting the Alps under tailwind conditions [(C) and (D)]. straight flight paths. (A to D) Tracks of migrating A. atropos in relation to Colors represent different individuals and are consistent between figures. Wind wind direction and speed (length of the arrows). The moths exhibit different layers are derived from the COSMO-1 model and represent conditions at 50 m strategies under different wind conditions, traveling due south through the Alps above ground level [(A) and (B)] and 300 m above ground level [(C) and (D)], the when primarily encountering headwinds [(A) and (B)] but traveling toward the estimated altitudes at which the moths were flying in the corresponding cases. The maintenance of consistently straight yet known of the capability of hawkmoths to Here, we provide evidence that large night- tracks and regulation of ground speed through- detect magnetic fields). At the landscape scale, flying insects actively select an orientation out the night under variable wind conditions we propose that the moths used topographical strategy in response to environmental con- strongly suggests that A. atropos has an inter- cues to visually navigate, because magnetic ditions, at least for some part of their mi- nal compass mechanism. Flight simulator cues are unlikely to be accurate enough to gratory journey. To maintain such straight studies have demonstrated that migrating maintain such straight trajectories. Overlaying trajectories over long periods of time, as seen Bogong moths (Agrotis infusa) use a combi- the straight tracks on a topographical map here, the moths must regularly update their nation of visual landmarks and Earth’s mag- (Fig. 1C) shows that the three orientation strat- position relative to whichever navigational netic field to navigate toward a goal (13). This has egies, and their directional clusters, are each cues they rely on. However, complete compen- yet to be demonstrated in free-flying migra- clearly aligned with a topographical feature that sation has not been previously documented in tory insects, but we predict that migrating would also result in avoiding the highest eleva- a long-range migratory insect and is generally hawkmoths, which have excellent nocturnal tions of the Alps (high-altitude passes running an unusual and very rare strategy in long- vision (24), use a similar suite of sensory mo- due south and southwest through the Alps and range migrants (25). Our results show that dalities to navigate over very large spatial a wide valley running west-southwest that would complex migratory strategies are not limited scales during migration (although nothing is enable circumventing the Alps altogether). to vertebrates. 766 12 AUGUST 2022 • VOL 377 ISSUE 6607 science.org SCIENCE

RESEARCH | REPORTS Fig. 3. Strategies of migrating hawkmoths in relation to winds. The flight behavior of migrating A. atropos in relation to winds encountered along the route was explored by analyzing distributions of the angle of deviation [b, the difference between the trajectory (T) and the wind direction (W)] for each segment of the trajectories of the individual moths shown in Fig. 1. (A) In unfavorable conditions such as headwinds and variable winds, moths had distributions of b with 95% CIs that overlapped 180° in all cases (table S1), indicating that they compensated for drift and maintained a southward track by selecting a flight heading (H) directly upwind. (B and C) Under favorable conditions (i.e., winds blowing toward the southwest), moths showed one of two strategies. As shown in (B), some moths had mean values of b around 45° and 95% CIs that did not overlap with 0° (table S1), indicating that they partially compensated for drift to migrate closer to south than the downwind flow would transport them. As shown in (C), moth 5 had a distribution of b that overlapped with 0° (table S1), indicating that it headed more-or-less straight downwind. On the circular plots, each point represents the value of b for all trajectory segments of each moth. The arrow indicates the overall mean value of b, and arrow length (r) indicates the degree of clustering around the mean. [Illustration credit: H. J. Williams] Fig. 4. Migrating hawkmoths modulate their airspeed in relation to wind represent the 25th and 75th percentiles [interquartile range (IQR)], and whiskers assistance. (A to G) Airspeed (m/s) of seven migrating A. atropos individuals [(A) to indicate data within ±1.5 times the IQR. The central bars represent median values. (G)] in relation to wind assistance (m/s) en route. Wind assistance was calculated Colors represent different individuals and are consistent between figures. Regression as the wind vector at the location of the moth in the direction of travel toward its next lines from linear models (LMs) are presented for significant relationships. LMs location, with positive values indicating tailwind and negative values indicating were performed for individuals with more than five data points. Significance headwind. Moths generally increased their airspeed under headwind conditions [(A) (P < 0.05) was based on likelihood-ratio tests: (A) moth 1, F = 651.91, P < 0.001; to (D)] and reduced their airspeed in more favorable tailwind conditions [(E) to (D) moth 11, F = 20.66, P < 0.001; (E) moth 4, F = 0.065, P = 0.807; (F) moth 6, (G)]. (H) Ground speed per segment for each of the seven individuals. Boxes F = 24.69, P = 0.008; and (G) moth 5, F = 0.465, P = 0.514. SCIENCE science.org 12 AUGUST 2022 ¥ VOL 377 ISSUE 6607 767

RESEARCH | REPORTS REFERENCES AND NOTES FLUID DYNAMICS 1. R. A. Holland, M. Wikelski, D. S. Wilcove, Science 313, 794–796 Dynamics of active liquid interfaces (2006). Raymond Adkins1†, Itamar Kolvin1*†, Zhihong You1†, Sven Witthaus1, 2. G. Hu et al., Science 354, 1584–1587 (2016). M. Cristina Marchetti1,2*, Zvonimir Dogic1,2* 3. K. R. Wotton et al., Curr. Biol. 29, 2167–2173.e5 (2019). 4. J. W. Chapman, D. R. Reynolds, K. Wilson, Ecol. Lett. 18, Controlling interfaces of phase-separating fluid mixtures is key to the creation of diverse functional soft materials. Traditionally, this is accomplished with surface-modifying chemical 287–302 (2015). agents. Using experiment and theory, we studied how mechanical activity shapes soft interfaces 5. D. A. Satterfield, T. S. Sillett, J. W. Chapman, S. Altizer, that separate an active and a passive fluid. Chaotic flows in the active fluid give rise to giant interfacial fluctuations and noninertial propagating active waves. At high activities, stresses P. P. Marra, Front. Ecol. Environ. 18, 335–344 (2020). disrupt interface continuity and drive droplet generation, producing an emulsion-like active 6. C. B. Williams, Insect Migration (Collins, 1958). state composed of finite-sized droplets. When in contact with a solid boundary, active 7. C. B. Williams, Proc. R. Entomol. Soc. London 13, 70–84 (1949). interfaces exhibit nonequilibrium wetting transitions, in which the fluid climbs the wall 8. C. G. Johnson, Migration and Dispersal of Insects by Flight against gravity. These results demonstrate the promise of mechanically driven interfaces for creating a new class of soft active matter. (Methuen, 1969). 9. L. R. Taylor, R. A. French, E. D. M. Macaulay, J. Anim. Ecol. 42, L iquid-liquid phase separation (LLPS) is a How cytoskeletal active stresses couple to self- ubiquitous phase transition, with exam- organization of membraneless organelles re- 751–760 (1973). ples abounding throughout material sci- mains an open question. Studies of simplified 10. J. W. Chapman et al., Science 327, 682–685 (2010). ence, biology, and everyday life (1, 2). systems can shed light on these phenomena. 11. J. W. Chapman et al., Curr. Biol. 18, R908–R909 (2008). Immiscible liquid phases are separated Relatedly, active wetting plays a potential role 12. J. W. Chapman et al., J. Anim. Ecol. 85, 115–124 (2016). by sharp but deformable interfaces that strong- in the development and shaping of tissues (18). 13. D. Dreyer et al., Curr. Biol. 28, 2160–2166.e5 (2018). ly couple to flows and the input of mechanical 14. D. Dreyer et al., J. Exp. Biol. 221, jeb179218 (2018). energy. For example, gentle shaking of an oil- To explore how activity modifies soft inter- 15. W. D. Kissling, D. E. Pattemore, M. Hagen, Biol. Rev. Camb. water mixture induces gravity-capillary inter- faces, we combined poly(ethylene glycol) (PEG) facial waves, whereas more vigorous perturba- and polysaccharide dextran with stabilized Philos. Soc. 89, 511–530 (2014). tions break up the entire interface, reinitializing microtubule filaments and clusters of kinesin 16. A. R. Pittaway, The Hawkmoths of the Western Palaearctic the phase separation (3–6). Active matter pro- molecular motors. Above a critical polymer vides an alternative method of continuously concentration, the passive PEG-dextran mix- (Harley Books, 1993). stirring a fluid (7, 8). In such systems, me- ture phase separated (19). Microtubules and 17. P. Howse, Bee Tiger: The DeathÕs Head Hawk-Moth Through the chanical energy, inputted locally through the kinesin clusters exclusively partitioned into motion of microscopic constituents, cascades the dextran phase, in which depletion forces Looking-Glass (Brambleby Books, 2021). upward to generate large-scale turbulent-like promoted microtubule bundling (Fig. 1, A to 18. G. F. McCracken et al., R. Soc. Open Sci. 3, 160398 (2016). dynamics (9–11). We studied how active stresses C). Streptavidin-bound kinesin clusters (KSA) 19. J. W. Chapman et al., Curr. Biol. 18, 514–518 (2008). and associated flows perturb soft interfaces and stepped along adjacent microtubules within 20. B. Gao et al., Proc. Biol. Sci. 287, 20200406 (2020). LLPS. Using experiment and theory, we iden- a bundle, driving interfilament sliding. The 21. Materials and methods are available as supplementary tified universal features of active-LLPS, includ- kinesin-powered bundle extensions contin- ing giant interfacial fluctuations, traveling uously reconfigured the filamentous network, materials. interfacial waves, activity-arrested phase sepa- generating large-scale turbulent-like flows, 22. R. B. Srygley, R. Dudley, Integr. Comp. Biol. 48, 119–133 (2008). ration, and activity-induced wetting transitions. similar to those previously studied (Fig. 1D) 23. K. J. Leitch, F. V. Ponce, W. B. Dickson, F. van Breugel, These results demonstrate how active matter (9). The PEG-dextran interfaces were suscepti- drives liquid interfaces to configurations that ble to large deformations by active stresses gen- M. H. Dickinson, Proc. Natl. Acad. Sci. U.S.A. 118, e2013342118 are not accessible in equilibrium. In turn, active erated within the dextran phase because of their (2021). interfaces are elastic probes that provide in- ultralow interfacial tension (<1 mN/m) (19). 24. A. L. Stöckl, D. C. O’Carroll, E. J. Warrant, Curr. Biol. 26, sight into the forces that drive active fluids— 821–826 (2016). for example, by allowing for the measurement We first visualized the phase separation dy- 25. J. W. Chapman et al., Curr. Biol. 21, R861–R870 (2011). of the active stresses. namics of active LLPS in ~30-mm-thick hori- 26. M. H. M. Menz et al., Data from: Individual tracking reveals zontal microscopy chambers. In such samples, long-distance flight-path control in a nocturnally migrating The active liquid interfaces we studied belong PEG-dextran interfaces had a nearly flat vertical moth. Movebank Data Repository (2022); https://doi.org/10. to a wider class of activity-driven boundaries profile (fig. S1). The quasi–two dimensional 5441/001/1.f4r24r5r. that includes lipid bilayers, colloidal chiral (2D) nature of the system was supported by a fluids, and interfaces between motile and im- nearly constant area fraction of the PEG-rich ACKNOWLEDGMENTS motile bacteria in a swarm (12–16). From a domains (fig. S2). In a passive system with biology perspective, LLPS has emerged as a microtubules but no kinesin motors, the drop- Special thanks to C. Ziegler for assistance in the field and providing ubiquitous organizational principle (2, 17). lets coalesced slowly (Fig. 1E and movie S1). images; the Swiss Ornithological Institute for providing access The addition of motors altered the coarsening to the field site at Col de Bretolet; M. Thoma for discussion 1Department of Physics, University of California at Santa kinetics. At intermediate KSA concentrations, regarding the migration of A. atropos; D. Dreyer and R. Massy for Barbara, Santa Barbara, CA 93106, USA. 2Graduate program active flows powered droplet motility, which discussion of circular statistics; I. Couzin for useful discussion of in Biomolecular Science and Engineering, University of increased the probability of droplets encoun- the project; and D. Dechmann, K. Safi, and two anonymous California at Santa Barbara, Santa Barbara, CA 93106, USA. tering each other and coalescing, thus speed- reviewers for constructive comments on the manuscript. Funding: *Corresponding author. Email: [email protected] (Z.D.); ing up coarsening dynamics (Fig. 1E and This project has received funding from the European Union’s [email protected] (M.C.M.); [email protected] (I.K.) movie S2). Higher KSA concentrations accel- Horizon 2020 research and innovation program under Marie †These authors contributed equally to this work. erated buildup of interfacial fluctuations, Skłodowska-Curie grant agreement no. 795568 to M.H.M.M. and the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) under Germany’s Excellence Strategy – EXC 2117 – 422037984 to M.W. Rothamsted Research receives grant-aided support from the Biotechnology and Biological Sciences Research Council of the UK. Author contributions: Conceptualization: M.H.M.M., M.W.; Data collection: M.H.M.M., H.-M.B.-S., M.W.; Data analysis: M.H.M.M., M.S.; Writing – original draft: M.H.M.M., J.W.C., M.W., H.J.W.; Writing – review and editing: M.H.M.M., H.J.W., H.-M.B.-S., D.R.R., M.S., J.W.C., M.W. Competing interests: The authors declare no competing interests. Data and materials availability: The datasets generated and analyzed during the current study are available in the Movebank Data Repository (26). Custom code for extraction of wind variables is available in data S1. License information: Copyright © 2022 the authors, some rights reserved; exclusive licensee American Association for the Advancement of Science. No claim to original US government works. https://www.science.org/about/science-licenses-journal-article-reuse SUPPLEMENTARY MATERIALS science.org/doi/10.1126/science.abn1663 Materials and Methods Figs. S1 and S2 Tables S1 and S2 References (27–41) Data S1 MDAR Reproducibility Checklist Submitted 22 November 2021; accepted 14 July 2022 10.1126/science.abn1663 768 12 AUGUST 2022 • VOL 377 ISSUE 6607 science.org SCIENCE

RESEARCH | REPORTS leading to an entirely different dynamical t as a function of the arc-length distance s along librium interfaces whose roughness is com- state in which droplets incessantly fused and parable with those measured at the lowest fissioned with each other (Fig. 1, E and G, and the interface (Fig. 2A). Interfacial fluctuations activities. movie S3). were described by time-averaged power spectra The dynamics of activity-driven interfacial To quantify the influence of activity on the fluctuations exhibited nontrivial spatiotem- PEG-dextran phase separation, we measured SðkÞ ¼  k2  with qk ¼ ∫dsqðs; tÞeÀiks. Be- poral correlations. To gather sufficient statis- the equal-time two-point correlation function kqk t, tics, we imaged ~10-mm-long active interfaces CðDr; tÞ ¼ hIðr þ Dr; tÞIðr; tÞir, where I = 1 in over a 2-hour interval. Space-time maps of local the dextran phase and –1 otherwise (fig. S3). cause of equipartition of thermal energy among interface height h(x, t) exhibited diagonally Spatial correlations decayed over a length scale streaked crests and troughs that were sug- x, defined by C(x) = 0.5, which is comparable Fourier modes, the spectrÀum of eÁquilibrium gestive of propagating waves (Fig. 2C). These with the average droplet size (fig. S4). For pas- interfaces is SðkÞ e Tk2= k2 þ k2c , where T translational modes were also evident in time- sive samples, x increased slowly in time (Fig. denotes tempperffiffiaffiffitffiffiuffiffiffirffiffieffi . The capillary wave lapse movies (movie S6). To characterize these 1F). Enhanced coarsening at intermediate KSA number kc ¼ Drg=g sets a crossover from modes, we measured the dynamic structure concentration was reflected by a much faster a gravity-dominated regime at large scales factor (DSF) of the interface height Dðk; wÞ ¼ initial growth of x than the passive case. At S(k) ~ k2 to a plateau at small scales, where ∫dxdteikxþiwthhðx′; t′Þhðx′ þ x; t′ þ tÞix′;t′ (Fig. high motor concentration, x peaked at ~1 hour surface tension attenuates fluctuations. Active 2D). Over a finite range of wave numbers, the and subsequently decayed to a finite plateau, DSF exhibited peaks at finite frequencies wp, xsteady. In parallel, average interface curvature interfacial fluctuations were markedly differ- confirming the presence of propagating modes k monotonically grew, reaching a sufficiently (Fig. 2E). Increased KSA concentration resulted large value to cause droplet fission (fig. S5). ent. Active spectra S(k) increased for small in higher wp for the same wave numbers; thus, The steady-state length scale xsteady was main- wave numbers (Fig. 2B and figs. S7 and S8). activity controlled the phase velocity (Fig. 2F). tained by the balance of droplet fission and After reaching a maximum for km ~ 30 mm–1, fusion events, in which xsteady was compara- it decayed as S(k) ~ k–3, instead of plateauing The giant nonequilibrium fluctuations and ble with the inverse of the average interface as in equilibrium. Whereas the shape of S(k) propagating wave modes result from the in- curvature ksteady. Concomitantly with the pla- remained the same for all KSA concentrations, teraction of active flows in the bulk dextran teauing of x, active flow speed became constant (Fig. 1F). These results demonstrate activity- the root mean square tangent angle qrms in- suppressed coarsening dynamics, which created creased linearly with activity (Fig. 2B, inset). an emulsion-like state in which finite-sized droplets continuously merge, break apart, and Using the crossover at km as a determinant exchange their content (Fig. 1G and movies S3 of the fluctuation amplitude, it would take an and S4). The volume fractions of the active and effective temperature of ~1011K to achieve equi- passive phases were nearly equal (fig. S2). Low volume fraction of active fluid generated sim- Fig. 1. Active LLPS. ilar steady states. Finite-sized domains are (A) Coexisting reminiscent of theoretical prediction in motility- PEG-rich (dark) and induced phase separation of isotropic active dextran-rich (cyan) particles (20). However, in contrast to theory, domains. (B) Labeled the active fluid in our experiments is aniso- microtubules (red) are tropic and perturbs an underlying equilibrium dissolved in the dextran phase separation. phase. Scale bar, 75 mm. (C) Microscopic- To gain insight into how active stresses drive scale depiction of phase interfacial fluctuations, we formed a macro- separation. Minority scopic interface through gravity-induced bulk PEG polymers (gray) in phase separation (Fig. 2A and fig. S6). In equi- the dextran-rich phase librium, molecular motion works against the induce microtubule density difference Dr and interfacial tension g bundling. (D) Kinesin to roughen the liquid-liquid interface. Typical clusters drive interfila- disturbances of PEG-dextran interfaces, bereft ment sliding. (E) Time of activity, are ~100 nm in amplitude, resulting evolution of the active in boundaries that appear flat when viewed LLPS at three KSA with our imaging setup (Fig. 2A and fig. S7). concentrations. Scale When driven out of equilibrium, however, bar, 350 mm. (F) (Top) interfaces exhibited giant undulations that Correlation length were visible with the naked eye (movie S5). evolution x(t) for As motor concentration increased, interfaces three KSA concentra- became multivalued with frequent overhangs, tions. For 230 nM KSA, indicating that active stresses directly control x plateaus at long interface configurations (Fig. 2A and movie S6). time (yellow highlight). (Bottom) Root mean The interplay of activity and capillarity is square velocity clarified by measuring the interfacial fluctu- of turbulent flows ation spectrum. To this end, local interface in the dextran tangent angles q(s, t) were sampled at a time phase at 230 nM KSA. (G) Fusion and fission of PEG droplets. Sample chamber thickness, 30 mm. Scale bar, 100 mm. SCIENCE science.org 12 AUGUST 2022 • VOL 377 ISSUE 6607 769