{'fraction_non_alphanumeric': 0.03367432378317376, 'fraction_numerical': 0.016382103462084532, 'mean_word_length': 5.040017590149517, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 3, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Improving AI System Awareness of Geoscience Knowledge: Symbiotic Integration of Physical Approaches and Deep Learning

Modeling dynamic geophysical phenomena is at the core of Earth and environmental studies. The geoscientific community relying mainly on physical representations may want to consider much deeper adoption of artificial intelligence (AI) instruments in the context of AI's global success and emergence of big Earth data. A new perspective of using hybrid physics‐AI approaches is a grand vision, but actualizing such approaches remains an open question in geoscience. This study develops a general approach to improving AI geoscientific awareness, wherein physical approaches such as temporal dynamic geoscientific models are included as special recurrent neural layers in a deep learning architecture. The illustrative case of runoff modeling across the conterminous United States demonstrates that the physics‐aware DL model has enhanced prediction accuracy, robust transferability, and good intelligence for inferring unobserved processes. This study represents a firm step toward realizing the vision of tackling Earth system challenges by physics‐AI integration.


Introduction
As a cornerstone of geosciences, physical approaches have achieved notable success in explaining and predicting the state changes in a geosystem (Bauer et al., 2015;Eyring et al., 2019), whereas nowadays they have to confront the development of artificial intelligence (AI), especially the recent development of deep learning (DL). In the context of the emergence of big Earth data (Yang et al., 2016), DL evolved into a budding tool for making scientific predictions independent of physical principles and its application is receiving a significant boost in various geoscientific domains (Brodrick et al., 2019;Hulbert et al., 2019;Rasp et al., 2018;Shen, 2018). Using AI techniques in geosciences has a long history. For example, the domain of hydroinformatics, formulated by Abbott (1991) 30 years ago, is defined as a union of computational hydraulics and AI (Solomatine & Ostfeld, 2008). However, the mainstream geoscientific community remains careful in adopting AI approaches, largely because of the assumed black-box nature of an AI model: it offers few mechanistic explanations beyond its fitting capability (Ebert-Uphoff et al., 2019; . Although some scientists have attempted to explain black-box models (Samek & Müller, 2019), Rudin (2019) recognizes that explaining black-box models rather than creating interpretable models in the first place is likely to perpetuate bad practice because the intentional explanations may be unreliable and even misleading. In addition to the reluctance in the community, data limitation is another obstacle (Karpatne et al., 2019).
Considering the respective merits of physical approaches and AI models, the synergy of the two paradigms has recently been envisioned as an attractive research topic in the geoscience community (Gil et al., 2019). Reichstein et al. (2019) summarized feasible strategies of physics-AI synergy in geoscience, some of which have been actively explored, such as using AI to build a surrogate model (Dwelle et al., 2019;Mo et al., 2019) or to learn and correct the mismatch between physical models and observations (Solomatine & Shrestha, 2009;. In comparison, the strategy of really hybrid modeling, which focuses on a more physically realistic neural network (NN) through adding one or several physical layers, appears to be much less investigated . The hybrid modeling approach more closely meets the prospect of improving AI systems' geoscientific awareness, since the entire sequence of tasks is undertaken with a single and unified AI architecture ( Figure S1 in the supporting information). Some studies in geoscience have attempted this strategy by introducing physical constraints into the loss function in the DL as a penalty term (Karpatne et al., 2017;Zhao et al., 2019) or wrapping analytical solutions of physical representations in DL as NN layers (de Bezenac et al., 2019;Wang et al., 2020). The existing approaches require physical knowledge to be expressed in a closed form. As the Earth is a highly dynamic system (Ghil, 2019), physical principles in geoscience are more often represented by differential equations (DEs) that cannot be analytically solved to achieve a closed-form solution. Consequently, how to design a network architecture for implanting geosystem dynamics when analytical solutions are not available remains an open question.
This study proposes an applicable solution to this problem by innovating a new NN architecture named physical process-wrapped recurrent neural network (denoted as P-RNN), which is designed to incorporate nonanalytically solvable ordinary differential equations (ODEs) into DL models. ODEs are extensively employed by geoscientists for studying geosystem dynamics such as seismic signal analysis (Peng et al., 2014) and global climate modeling (Kaper & Engler, 2013) since the nineteenth century (Ghil, 2019). In this study, the new architecture was tested in the field of hydrology, using runoff modeling across the conterminous United States (CONUS) as an illustrative case. Hydrology is a cornerstone discipline in the geosciences (Fatichi et al., 2016), and catchment runoff modeling is the core task of hydrology (Beven, 2012). A catchment is a typical dynamic geosystem where water input (e.g., rainfall), output (e.g., runoff), and storage (e.g., soil moisture and snowpack) evolve over time. In this illustrative case, a conceptual hydrologic model (EXP-HYDRO; Patil & Stieglitz, 2014) is included as the P-RNN layer in a DL architecture, leading to hydrology-aware DL models. Enhanced simulation accuracy, robust cross-catchment transferability, and good intelligence for inferring unobserved processes were observed for the hydrology-aware DL models, highlighting the symbiotic integration between DL and the physical approach: the data-driven components fill gaps in the physical knowledge underlying the hydrologic model and the knowledge imparts physical awareness to the DL model. Overall, this study demonstrates that AI can garner physical knowledge the way humans do, if taught appropriately, and a symbiotic integration between DL and physical approaches is feasible, beneficial, and promising in advancing geoscience.

Generic Architecture
A physics-aware DL architecture comprising two data pipelines, each with several NN layers, is proposed (see Figure 1a). In the main pipeline, between the input (i.e., dynamic forcings) and output (i.e., dynamic responses), a P-RNN layer (to be elaborated in section 2.2) is constructed to wrap the geoscientific model (block a). The outputs of the P-RNN layer enter a series of common NN layers (block b), which can be fully connected layers, convolution layers, or recurrent layers (Goodfellow et al., 2016). In this design, the P-RNN layer endows the network with physical interpretability and consistency, and the common NN layers address processes unrepresented by the P-RNN layer (when the output of the P-RNN layer is different from the final output) or correct mismatch(when the output of the P-RNN layer is the same as the final output).
The main pipeline contains two types of parameters: the architecture parameters θ N (i.e., weights and biases) of the common NN layers and physically meaningful parameters θ P of the P-RNN layer. To improve the transferability of the physics-aware DL across different regions, the architecture involves a parameterization pipeline for mapping region-dependent attributes (e.g., soil properties and topography) onto θ N and θ P using additional networks (blocks c and d). Once the architecture parameters of the parametrization pipeline itself are obtained, θ N and θ P can be determined on the interior. The unique design of a separate parameterization pipeline makes the architecture flexible to produce either global/regional models (the pipeline enabled) or local models (the pipeline disabled) and ensures that the architecture can infer spatially varying θ N and θ P when addressing cross-region learning.

The P-RNN Layer
The proposed P-RNN layer is designed to encode geosystem dynamics as an NN architecture. A dynamic geosystem is often characterized by the state-space representation, which describes the system response for certain inputs and can be modeled by a combination of ODEs (state equations) and output equations in the following form: where s(t) ∈ R n , x(t) ∈ R p , and y(t) ∈ R q , respectively, denote state vectors, input vectors, and output vectors with respect to time t, F : R n × R p → R n , G : R n × R p → R q , and θ 〈F〉 and θ 〈G〉 are the parameter sets.
Analogous to the ordinary recurrent neural network (RNN) architecture (Rumelhart et al., 1986), the backbone of the P-RNN layer is formed by recurrent cells that can provide memories of the past sequence (see Figure S2). Within recurrent cells in the P-RNN architecture, the connections between neurons (inputs x, states s, and outputs y) are specified with the state-space representation (Equation 1) in an explicit discrete form and the architecture parameters (i.e., weights and biases) of an ordinary RNN are replaced by parameters with physical meanings (i.e., θ P ). Table S1 in the supporting information summarizes the pseudocode of the P-RNN layer implemented in this study. Niu et al. (2019) demonstrated the connections between network architectures of the RNN family and numerical methods of ODEs, which theoretically supports the use of the P-RNN to tackle problems involving system dynamics.

Implementation for Runoff Modeling
In the case study, the EXP-HYDRO model (Patil & Stieglitz, 2014) is wrapped with the P-RNN layer. It is a conceptual, spatially lumped hydrologic model underpinned by two discrete state-space representations ( Figure 1b). Although the model is not process-based (or physically based) per the rigorous definition in Fatichi et al. (2016), it adheres strictly to the law of mass conservation and has been demonstrated competent in large-sample catchments (Patil & Stieglitz, 2014). The snow accumulation bucket (S 0 ) and the storage in the catchment bucket (S 1 ) are the model state variables. Daily precipitation P, daily temperature T, and day length L day are the model inputs, and the daily streamflow at the catchment outlet Q is the output of interest. The model involves six physically meaningful parameters (i.e., θ P ) that control the hydrologic behaviors (see Figure 1b). Detailed information on the model is provided in Text S1 in the supporting information.
In this case, the two pipelines in the generic architecture ( Figure 1a) are respectively for runoff modeling and its parameterization. The inputs of the main pipeline are the meteorological forcing variables including P, T, and L day , as required by the P-RNN layer, as well as the shortwave downward radiation SRad and vapor pressure VP. The EXP-HYDRO model-wrapped P-RNN layer provides a preliminary runoff estimation Q * , which is fed into a two-layer common NN block together with the five input variables (see Figure S3). In this study, the 1-D convolution layers (Conv1D; Fukushima & Miyake, 1982) are utilized in the common NN block for considering the lagged influence of predictors in 10 days on current hydrologic responses (see Text S2 for details). The Conv1D layer is capable of handling the lagged effect by a one-direction convolution operation and has been employed for data-based hydrologic modeling in several studies (Feng et al., 2019). Through the Conv1D layers, the approximation errors associated with the physical approach are corrected, and the final runoff Q is derived. The model following the main pipeline scheme that hybridizes physical approaches (represented by the P-RNN layer) and data-driven components (i.e., Conv1D layers) is hereafter referred to as the "hybrid DL model." Furthermore, the parameterization pipeline is used to provide the catchment awareness for the main pipeline, in which two blocks of fully connected layers offer θ P for the P-RNN layer and θ N for the Conv1D layers. As a result of the parameterization pipeline, θ N and θ P can vary among different catchments with the physiographic attributes. More information regarding the network structures (e.g., number of layers and neurons) and training hyperparameters (e.g., training epoch and learning rate) of the model is elaborated in Text S2.
The data used in this study were obtained from the freely available CAMELS dataset (Addor et al., 2017;Newman et al., 2015), which consists of the area-averaged daily hydrometeorological time series and catchment attributes for 671 catchments across the CONUS with limited human disturbance. In this study, 569 catchments where daily hydrologic observations are complete and continuously recorded from 1 October 1980 to 30 September 2010 were considered (see Figure S4).

Synergy of Physics and AI
This study first built separate hybrid DL models for individual catchments, in which only the main pipeline in Figure 1a was activated, as each model only involves a single catchment where the spatial heterogeneity of parameters is beyond consideration. For each individual catchment, the hybrid DL model was compared against three other models. The first two models are variants (two partial models). The first variant is referred to as the "physical NN model," in which the main pipeline contains the "physical" part only (i.e., the P-RNN layer); the second variant is referred to as the "common NN model," in which the main pipeline only contains the "data-driven" part (i.e., the Conv1D layers). The third model is the long short-term memory (LSTM) network hydrologic model developed by Kratzert et al. (2018) for reference purposes. Both the common NN model and the LSTM model are purely data-driven, whereas the physical NN model is essentially a numeric variant of the EXP-HYDRO model. The four models were separately calibrated at 569 individual catchments for daily streamflow from 1 October 1980 to 30 September 2000. Figure 2a presents the spatial distributions of the Nash-Sutcliffe efficiency (NSE) values, which measure the overall simulation performance (Nash & Sutcliffe, 1970), from the four models during the independent evaluation period (1 October 2000 to 30 September 2010). Visually, the hybrid DL model outperforms the other three in terms of the NSE in general. For all the models, the catchments with low NSE values are mainly located in the central arid regions of the United States in which runoff is largely controlled by short-duration, high-intensity precipitation events and infiltration-excess overland flow dominates runoff generation (Berghuijs et al., 2016). The prototypical model of the P-RNN layer operates with a mechanism of saturation-excess overland flow and is mainly applicable in the Northeast and the Pacific Northwest as well as forested mountain areas in the United States, which is indicated by the worst performance of the physical NN model outside these regions. It is also difficult for purely data-driven models such as the common NN model and the LSTM model to satisfactorily reproduce flashy hydrographs in such regions due to the infrequent storms and flood records. Such a spatial pattern of model performance was also revealed in previous studies (Essou et al., 2016;Herman et al., 2013;Kratzert et al., 2018;Newman et al., 2017;Sun et al., 2014) with similar explanations (McCabe & Wolock, 2011).
We further compared the models by scrutinizing the hydrographs they reproduce for the individual catchment. The hybrid DL model shows superior performance in reproducing flow peaks and capturing the overall pattern. The two purely data-driven models show unfavorable performance for flow peaks, and the physical NN model demonstrates mediocre overall performance. Figure S5 shows a detailed comparison of the hydrographs predicted by the four models at a typical catchment. Similar results were observed for the majority of other catchments, as indicated by Figure 2b, which illustrates the distributions of three goodness-of-fit metrics, including the NSE, the absolute value of peak flow bias (PFAB; Yilmaz et al., 2008), and mNSE (emphasizing the performance for baseflow; Legates & McCabe, 1999) evaluated for the four models. Table S2 compares the NSE performance of the four AI models in this study with those of previous models (Kratzert et al., 2018;Newman et al., 2017;Patil & Stieglitz, 2014;Xia et al., 2012), most of which are based on the same dataset and have an overlapping evaluation period with that in this study. Among the models compared, the hybrid DL model exhibits the highest overall performance. Table S3 explores the correlation between the hybrid DL model and its two partial models in terms of the NSE performance within the 569 catchments. NSE ≥ 0.55 was considered the threshold for good performance (Knoben et al., 2019;Newman et al., 2015). If at least one of the partial models achieves good performance (see rows 1 to 3 and 5 to 7 of the table), the hybrid DL model performs satisfactorily in most cases, which implies that the predictive power of the partial models can be maintained by the hybrid model. More interestingly, even when both partial models do not meet the performance criterion, there is still a large chance (see row 4 of the table) that the hybrid model reaches a satisfactory NSE, which reminds us of the old saying "two heads are better than one." Further investigation reveals that the hybrid model can adequately maintain water balance when the prototypical model encoded in the P-RNN layer has a reliable performance in the target catchment (see Figure S6).

Robust Transferability of the Educated AI
Hydrologic predictions in ungauged basins are a great challenge in hydrology (Hrachowitz et al., 2013) and require better comprehension of the links between the hydrological function and physical properties of a catchment (Wagener et al., 2007). A typical way to address predictions in ungauged basins is to build a regional hydrologic model with site-specific parameterization (Oudin et al., 2008;Razavi & Coulibaly, 2013). In this study, a regional DL model is built following the new DL architecture that contains both the main pipeline and parameterization pipeline and is hereafter referred to as the "regional hybrid model." As mentioned above, runoff modeling based on the CAMELS dataset is generally unsuccessful in the arid central part of the United States. Figure S7 suggests the large possibility that the hybrid DL model performs satisfactorily for an aridity index ≤1.2 or the frequency of dry days ≤280 days/year. Thus, 450 catchments (see Figure S8) were retained to evaluate the model's transferability via fivefold cross validation. Specifically, the 450 catchments were randomly and equally sampled into five groups. In each of the five runs, the regional hybrid model, which employs 27 available physiographic attributes (summarized in Table S4) as the inputs for parameterization, was trained during 1 October 1980 to 30 September 2010 in four groups of catchments (as gauged) and validated on the rest of the groups (as ungauged) during the same period. For comparison, two additional regional models were built and trained in the same way. One is the same as the regional hybrid model except that the parameterization pipeline is not activated, and therefore, θ P and θ N are constant across catchments. This model is hereafter referred to as the "pseudo-regional hybrid model." The other is the regional LSTM model by Kratzert et al. (2019), which uses the direct concatenation of catchment attributes and meteorological variables as inputs. The pseudo-regional hybrid model ignores the catchment attributes, and the regional hybrid model and regional LSTM model represent two distinct strategies of representing catchment attributes. Figure 2c illustrates the box plots of the NSE values achieved by the three regional models in the cross validations. Compared with the NSE values from the pseudo-regional hybrid model, incorporating catchment attributes in the regional hybrid model significantly improves the performance, which confirms that the variations in physiographic attributes between catchments are critical to distinguishing hydrologic behaviors. Although the regional LSTM model exhibits slightly better performance than the regional hybrid model with respect to the median NSE, its predictive ability has significant fluctuations and the model showed catastrophically bad predictions (NSE < 0) in 23 catchments. The mean NSE achieved by the regional hybrid model is higher, the standard deviation is considerably lower (an important indication of a robust model), and fewer catastrophic predictions occur (see Figure S9), suggesting a much more robust transfer-learning ability of the regional hybrid model.
LSTM models are currently the popular machine learning approach for runoff modeling. The regional LSTM model in this study infuses the catchment attributes and meteorological forcing into one neuron perceptron (these values are directly weighted and summed) and achieved a slightly higher median NSE than the regional hybrid model, whereas the interpretable structure of the latter enhances the confidence of predictions in unfamiliar regions. As noted by Hrachowitz et al. (2013), linking catchment form (i.e., the attributes) to hydrological function (i.e., the way a catchment responds to input, including model parameters) is important for adequately addressing the challenge of ungauged basins. From this perspective, it is worth giving priority to the model's inherent interpretability even if it allows for some decrease in accuracy. Moreover, by dissecting the parameterization pipeline in a calibrated model, one can further explore spatial patterns of the inferred parameters, which offer insights into catchment processes and performance of the hydrologic model (e.g., Beck et al., 2013;Foks et al., 2019;Santhi et al., 2008). In this sense, the regional hybrid model can also be viewed as a streamlined parameterization module for the original hydrologic model, which is highly efficient in cross-region applications. Please refer to Text S3 and Figure S10 for more details.

Can the Educated AI Reason?
The previous sections show that the symbiotic integration between DL and physical approaches can realize better and more robust predictive power. However, the predictive power mainly reflects the induction capability of the educated AI rather than its deduction capability. The hydrologic model learned by the AI system describes runoff generation as well as the dynamics of snow accumulation and melting, and the P-RNN layer can output the water storage in snowpack (denoted as S 0 ) as an intermediate state variable. Therefore,

Geophysical Research Letters
it is of great interest to determine whether the AI system trained with flow observations only can appropriately infer snow pattern and dynamics. An affirmative answer would serve as evidence for the reasoning capability of educated AI. Figure 3a illustrates the average S 0 from 1 October 2000 to 30 September 2010 in all 450 catchments based on the regional hybrid model simulations. Intuitively and without exception, the catchments with noticeable snow accumulation have high altitudes or high latitudes. To verify the unsupervised learned S 0 , daily snow water equivalent (SWE) data over the CONUS (4 km × 4 km spatial resolution) were obtained from the NASA National Snow and Ice Data Center (Zeng et al., 2018). Figure 3b presents the spatial distribution of the average SWE in the same period, and the spatial pattern coincides with that in Figure 3a. More Figure 3. Comparison of the water storage in snowpack (S 0 ) inferred by the regional hybrid model with the snow water equivalent (SWE) data. (a) Ten-year average of daily S 0 in 450 catchments obtained from the P-RNN layer in the regional hybrid model. The highlighted catchments with letters were selected for the analyses in Figure S11. surprisingly, as Figure 3c shows, the snow dynamics inferred by the educated AI perfectly match those of the independent SWE dataset (with a correlation coefficient of 0.99). Similar comparisons were also made for individual catchments. Figure S11 demonstrates an excellent match for representative catchments (highlighted in Figure 3a). These promising results confirm that the DL model with physical awareness is capable of not merely undertaking the intentional task of data fitting, but more importantly, reasonable deduction of processes without direct observations.

Toward Smarter AI for Geoscience
This study seamlessly merged physical knowledge into DL architecture by creating a P-RNN layer to solve ODEs in discrete form, a common way to describe geosystem dynamics. It is important to note that the architecture presented is not necessarily an unchangeable template, considering the flexibility and easy extensibility of DL (e.g., modular building blocks). For example, the common NN layers noted in Figure 1a can be any type, such as the advanced ConvLSTM layers (Shi et al., 2015) that capture spatiotemporal correlations.
Purely data-driven NN layers can also be put on top of the physics-aware NN layers to provide inputs (e.g., parameterizations for clouds (Rasp et al., 2018) and the fluid motion field (de Bezenac et al., 2019)).
In addition to solving ODEs with NN layers to reproduce 1D system dynamics, encoding partial differential equation (PDE) solvers in DL architectures to simulate higher-dimensional system dynamics deserves further exploration in geoscience. In recent years, DL has gradually emerged as a powerful technique for solving PDEs in the field of applied mathematics (Sirignano & Spiliopoulos, 2018). For example, the utilization of NN layer architectures involving iterative convolution operations is a promising solution for the combination of DL and higher-dimensional system dynamics (Long et al., 2018).
With regard to future applications in runoff modeling, one of the next steps would be to encode additional runoff generation mechanisms and hydrologic processes (e.g., subsurface flow) in the P-RNN layer, such that a higher transferability of regional models can be achieved.

Conclusion
This study represents a firm step toward realizing the vision of Reichstein et al. (2019) of tackling Earth system challenges using hybrid physics-AI approaches. We demonstrate that symbiotic integration is critical to the success of such approaches, as the data-driven component fills gaps in physical knowledge (e.g., processes unrepresented by physical approaches and physical parameters that are difficult to measure) and the physical component endows the AI system with physical consistency and reasoning ability. As the illustrative case of runoff modeling shows, the DL architecture proposed in this study enables an AI system to directly learn elements of physical representations and make accurate deductions on unobserved phenomena. In this regard, the novel framework presents a promising strategy to effectively educate AI with available geoscientific knowledge. We envision a number of future studies to implement the developed framework in various geoscientific contexts and further expand it to enable the use of smarter AI systems to advance research in the geoscience field.