Výsledky a diskuse
Autoionizace byla zkoumána pomocí ab initio simulací RETIS, jak je popsáno v části Materiály a metody. Pro simulace RETIS jsme použili relativně jednoduchý geometrický parametr pořadí vzdáleností, λ, jak je znázorněno na obr. 1: Pokud se systém skládá pouze z druhů H2O, λ je největší vzdálenost kovalentní vazby O-H, a pokud systém obsahuje druhy OH- a H3O+, λ je brána jako nejkratší vzdálenost mezi atomy kyslíku v OH- a vodíku v H3O+. V dalším textu označujeme atom kyslíku použitý pro parametr uspořádání jako Oλ. Typ druhu (OH-, H2O nebo H3O+) byl určen tak, že každému vodíku byla přiřazena jedna vazba spojující jej s nejbližším kyslíkem. Všimněte si, že definice parametru uspořádání nevyžaduje prahovou hodnotu pro definování chemické vazby ani neomezuje parametr uspořádání na konkrétní molekuly vody po dobu trvání simulace. To znamená, že počítáme rychlost disociace jakékoli molekuly vody v systému namísto jediné cílené vazby O-H nebo molekuly vody.
Iniciační podmínky a lokální kolektivní proměnné. (A) Reaktivní (rλc,λr(ξ)) a nereaktivní (uλc,λr(ξ)) rozdělení pro ξ={w4,na} a λc=1,16 Å a λr=2,0 Å. Pro účely vizualizace jsou zobrazená rozdělení normalizována . Horní a pravá vložka zobrazují jednorozměrné projekce rozdělení. Je vidět jasné oddělení obou rozdělení podél souřadnice w4, což naznačuje, že reaktivní trajektorie jsou ve srovnání s nereaktivními trajektoriemi více stlačené. Kromě toho atom kyslíku použitý při výpočtu parametru uspořádání (Oλ) přijímá v průměru větší počet vodíkových vazeb v reaktivních trajektoriích ve srovnání s nereaktivními trajektoriemi. (B) Ilustrativní snímek z reaktivní trajektorie, kde je Oλ zobrazen modře. Oranžově jsou znázorněny čtyři okolní atomy kyslíku, které se používají pro výpočet tetraedrického parametru uspořádání q. Vodní vodič je zvýrazněn žlutou čarou (a šedými průhlednými kuličkami) a je vyznačen úhlový parametr qcos. Na tomto snímku je vodní drát stlačený, q vykazuje odchylku od tetraedrické struktury, qcos naznačuje, že v drátu jsou seřazeny tři atomy kyslíku a Oλ přijímá tři vodíkové vazby a jednu daruje (znázorněno zelenými čarami).
Pokud vezmeme v úvahu souřadnici q, zjistíme, že rλc,λr je posunuto směrem k nižším hodnotám q ve srovnání s uλc,λr, což naznačuje, že odklon od tetraedrického uspořádání kolem disociující vodní formy může také iniciovat událost. Toto zjištění je poněkud překvapivé, protože v některých jiných chemických reakcích ve vodné fázi byl zjištěn opačný efekt (31). Podobné závěry lze vyvodit i pro rozložení ξ=(w4,qcos). Zde se objevuje vrchol podél souřadnice qcos pro reaktivní rozložení bližší lineárnímu uspořádání molekul vody. Na obr. 4B ukazujeme reprezentativní snímek získaný na počátku (po 3 fs) reaktivní trajektorie. Celkově výsledky uvedené na obr. 3 sdělují, že stlačení vodního drátu (měřené pomocí w4) a hyperkoordinace (měřená pomocí na) nebo deformace (měřená pomocí q a qcos) jsou nezbytnými iniciačními podmínkami pro autoionizaci. Nejsou to však postačující podmínky, jak ukazují hodnoty TAλc,λr na obr. 3B: stále ještě 60 % trajektorií začínajících v ideálním rozsahu parametrů ξ nedokáže vytvořit souhlasný protonový skok.
Machine learning (ML) aplikovaný na data ze vzorkování cest (33, 34) je slibný přístup k nalezení důležitých kolektivních proměnných, které mohou být snadno přehlédnuty lidskou intuicí. Abychom tuto možnost prozkoumali, sestavili jsme ML modely pro předpovídání výsledku trajektorií vzhledem ke stavu vodního systému na počátku trajektorií. Zaměřili jsme se na stejný rozsah jako v analýze predikční síly a k předpovědi výsledku jsme použili stav systému, kdy je poprvé dosaženo λ>1,15 Å . Použili jsme několik technik ML, při nichž byl každý soubor lichých drah zahrnut do kalibrace a soubory sudých drah byly použity pro testovací soubor. Alternativní rozdělení, při kterém byla data v rámci každého souboru cest rovnoměrně rozdělena na dvě části, poskytlo podobné výsledky. Protože silně zkreslená rozdělení je navíc obtížné zpracovat pomocí ML, vynechali jsme dále převážení datových souborů statistickými váhami příslušných souborů cest. Techniky ML jsme však použili jako kvalitativní přístup k nalezení nových parametrů, které by mohly být testovány kvantitativně v rámci metody prediktivní síly (19).
Abychom se vyhnuli potenciálnímu riziku nadměrné interpretace, rozhodli jsme se navíc omezit složitost rozhodovacího procesu ML a při výpočtu TAλc,λr jsme stanovili maximální počet čtyř řádových parametrů. Vynikajících predikčních výkonů (>90 %) bylo dosaženo například při použití strojů založených na ensemble gradient-boosting (35, 36). Interpretace modelu je však problematická, protože se používá ansámbl 100-150 hlubokých rozhodovacích stromů (přidávaných postupně). Výkonnost se sice zlepšuje, ale zvyšuje se pravděpodobnost nadměrného přizpůsobení s náhodnými korelacemi. Omezili jsme se proto na rozhodovací modely založené na jednom stromu založené na klasifikačních a regresních rozhodovacích stromech (CART) (20). Omezení na čtyři řádové parametry funkce TAλc,λr vychází z podobných důvodů. Přidáním více parametrů získáme řidší matice reprezentující reaktivní/nereaktivní rozdělení, v důsledku čehož se numerická integrace pro výpočet překryvu mezi těmito rozděleními stává velmi citlivou na velikost binů a mohla by podhodnotit překryv v důsledku toho, že jsou biny prázdné kvůli nedostatečné statistice.
Uvažovali jsme 138 kolektivních proměnných sestávajících ze vzdáleností kyslík-kyslík; vzdáleností kyslík-vodík pro původně vázané molekuly vody; všech úhlů tvořených Oλ a jeho čtyřmi nejbližšími sousedy kyslíku; a Steinhardtovy parametry řádu 3, 4 a 6 (32) (podrobněji viz Materiály a metody). Kromě toho byly přidány již uvažované řádové parametry. Obr. 5A ukazuje výsledný rozhodovací strom. Pozoruhodné je, že ze všech vstupních parametrů je parametr w4 na vrcholu rozhodovacího stromu a zároveň nejdůležitější proměnnou, měřeno snížením klasifikační chyby připisované každé proměnné při každém rozdělení rozhodovacího stromu (20) (příloha SI, obr. S9). V rozhodovacím stromu se objevuje také tetraedrické uspořádání a počet akceptovaných vodíkových vazeb. Pro popis prvního efektu ML přístup upřednostnil Steinhardtův parametr pořadí q4 před podobným parametrem q, který jsme dříve používali my. Některé vzdálenosti, které se rovněž objevují v rozhodovacím stromu, jako d25, vzdálenost mezi Oλ a jeho 25. nejbližším kyslíkem, jsou pravděpodobně způsobeny náhodnými korelacemi způsobenými omezenou velikostí souboru dat. To lze ověřit kontrolou důležitosti této proměnné: d25 se neobjevuje mezi 20 nejdůležitějšími proměnnými (příloha SI, obr. S9) a ve skutečnosti jsou jiné podobné proměnné (např. d24) řazeny výše, i když s nízkou důležitostí. Důležitějším a intuitivně správným parametrem, který navrhuje ML přístup, je λ2, vzdálenost OH mezi kyslíkem nejblíže Oλ a jeho vodíkem s největší intramolekulární vazbou. Přepočet predikční schopnosti s použitím parametrů z ML stromu (obr. 5B) nepřinesl vyšší výkonnost než kombinace w4, q, na a qcos, ale měl by být chápán jako stejně dobrý s ohledem na statistické nejistoty.
Výsledky analýzy strojového učení. (A) Klasifikační a regresní strom pro předpověď výsledku iniciovaných trajektorií. Zde jsme uvažovali několik dalších kolektivních proměnných (popis v Materiálech a metodách), ale pro konstrukci stromu je nakonec potřeba jen malá podmnožina: w4, q4 , λ2 (délka natažené vodíkové vazby v molekule vody nejbližší druhu Oλ), di (vzdálenost Oλ od i-tého nejbližšího kyslíku) a d¯i (průměrná vzdálenost s ohledem na i nejbližších kyslíků). Zápis uzlů je vysvětlen pomocí samostatného uzlu v levém horním rohu. Tento strom předpovídá trajektorie jako reaktivní, tj. dosahující λ≥2, nebo nereaktivní na základě kolektivních proměnných získaných v rámci v trajektoriích, kdy λ je poprvé ≥1,15. Uzly předpovídající reaktivní trajektorie jsou podbarveny modře (třída 1), zatímco uzly předpovídající nereaktivní trajektorie jsou podbarveny zeleně (třída 0). Všimněte si, že procenta ve spodní části čtverců neodrážejí fyzikálně správné podíly, protože soubory drah nebyly převáženy pomocí svých statistických vah. Pravidla jsou textovou reprezentací procházení stromu; například pravidlo 5 (které předpovídá reaktivní trajektorie) lze vyjádřit jako w4≥7,6 a λ2≥1,1. Tato pravidla udávají různé iniciační podmínky a jsou uvedena v příloze SI, tabulce S1, pro spodní řadu uzlů. (B) Predikční síla a pravděpodobnost křížení v závislosti na λr pro λc=1,16 Å a různé kombinace kolektivních proměnných. Zde porovnáváme predikční sílu pomocí námi identifikovaných kolektivních proměnných s proměnnými, které byly analýzou strojového učení označeny jako důležité. (C) Reaktivní (rλc,λr(ξ)) a nereaktivní (uλc,λr(ξ)) rozdělení pro ξ={λ2,d¯2} a λc=1,16 Å a λr=2,0 Å. Pro účely vizualizace jsou zobrazená rozdělení normalizována. Horní a pravá vložka zobrazují jednorozměrné projekce rozdělení.
.