Author Markus ReichsteinAuthor Lazaro AlonsoAuthor Bernhard Ahrens

Soil respiration

using EasyHybrid
using CairoMakie

Create synthetic data

Exponential temperature relationship

using Random
Random.seed!(2314)

T = rand(500) .* 40 .- 10      # Random temperature
SM = rand(500) .* 0.8 .+ 0.1   # Random soil moisture
SM_fac = exp.(-8.0*(SM .- 0.6) .^ 2)
Resp0 = 1.1 .* SM_fac # Base respiration dependent on soil moisture
Resp = Resp0 .* exp.(0.07 .* T)
Resp_obs = Resp .+ randn(length(Resp)) .* 0.05 .* mean(Resp);  # Add some noise

df = DataFrame(; T, SM, SM_fac, Resp0, Resp, Resp_obs)
first(df, 5)

5×6 DataFrame

Row	T	SM	SM_fac	Resp0	Resp	Resp_obs
	Float64	Float64	Float64	Float64	Float64	Float64
1	25.3743	0.182296	0.247631	0.272394	1.60913	1.71876
2	-2.16137	0.798355	0.729965	0.802961	0.69022	0.733988
3	-0.854465	0.37379	0.66407	0.730477	0.688066	0.903229
4	14.2647	0.148261	0.195431	0.214974	0.5835	0.63434
5	5.41821	0.797535	0.731863	0.805049	1.17636	1.14584

Define the Process-based Model

function Expo_resp_model(;T, Resp0, k)
    Resp_obs = Resp0 .* exp.(k .* T)
    return (; Resp_obs, Resp0, k)
end;

about Expo_resp_model

Mechanistic model for soil respiration based on exponential temperature response

Mathematical Model: Resp = Resp0 * exp(k * T). Implements Arrhenius-type temperature dependence for biological processes (true connection?).

Arguments

• T : Air or soil temperature [°C]

• Resp0 : Basal respiration rate - which is equal to Resp at T=0°C [μmol CO₂ m⁻² s⁻¹]

• k : Temperature sensitivity parameter [°C⁻¹], typical range 0.05-0.15

Returns

• Resp_obs : Predicted soil respiration [μmol CO₂ m⁻² s⁻¹]

• Resp0 : Basal respiration (passed through)

• k : Temperature sensitivity (passed through)

Scientific Context

Models temperature dependence of soil CO₂ efflux for carbon cycle studies. Exponential relationship reflects enzyme kinetics in microbial decomposition.

Configure Model Parameters

parameters = (
    # name: (default, lower_bound, upper_bound) # Description
    k     = (0.01f0, 0.0f0, 0.2f0),  # Exponent
    Resp0 = (2.0f0,  0.0f0, 8.0f0),  # Basal respiration [μmol/m²/s]
);

Construct the Hybrid Model

targets = [:Resp_obs]
forcings = [:T]
predictors = (Resp0=[:SM],);

Define global parameters (none for this model, Q10 is fixed)

global_param_names = [:k]

1-element Vector{Symbol}:
 :k

hybrid_model = constructHybridModel(
    predictors,
    forcings,
    targets,
    Expo_resp_model,
    parameters,
    global_param_names,
    scale_nn_outputs=false, # TODO `true` also works with good lower and upper bounds
    hidden_layers = [16, 16],
    activation = sigmoid,
    input_batchnorm = true
)

Hybrid Model (Multi NN)
Neural Networks:
  Resp0:
    Chain(
        layer_1 = BatchNorm(1, affine=false, track_stats=true),
        layer_2 = Dense(1 => 16, σ),                  # 32 parameters
        layer_3 = Dense(16 => 16, σ),                 # 272 parameters
        layer_4 = Dense(16 => 1),                     # 17 parameters
    )         # Total: 321 parameters,
              #        plus 3 states.
Configuration:
  predictors:
    Resp0 = [:SM]
  forcing = [:T]
  targets = [:Resp_obs]
  mechanistic_model = Expo_resp_model
  neural_param_names = [:Resp0]
  global_param_names = [:k]
  fixed_param_names = Symbol[]
  scale_nn_outputs = false
  start_from_default = true
  config = (; hidden_layers = [16, 16], activation = σ, scale_nn_outputs = false, input_batchnorm = true, start_from_default = true,)

Parameters:
  Hybrid Parameters
    ┌───────┬─────────┬───────┬───────┐
    │       │ default │ lower │ upper │
    ├───────┼─────────┼───────┼───────┤
    │     k │    0.01 │   0.0 │   0.2 │
    │ Resp0 │     2.0 │   0.0 │   8.0 │
    └───────┴─────────┴───────┴───────┘

Train the Model

out =  train(hybrid_model, df, (:k,); nepochs=300, batchsize=64,
    opt=AdamW(0.01, (0.9, 0.999), 0.01), loss_types=[:mse, :nse],
    training_loss=:nseLoss, random_seed=123, yscale = identity,
    monitor_names=[:Resp0, :k],
    show_progress=false,
    hybrid_name="expo_response"
    );

[ Info: Training data type: AxisKeys.KeyedArray{Float32, 2, NamedDims.NamedDimsArray{(:variable, :batch_size), Float32, 2, SubArray{Float32, 2, LinearAlgebra.Transpose{Float32, Matrix{Float32}}, Tuple{Vector{Int64}, UnitRange{Int64}}, false}}, Tuple{SubArray{Symbol, 1, Vector{Symbol}, Tuple{Vector{Int64}}, false}, UnitRange{Int64}}}
[ Info: Check the saved output (.png, .mp4, .jld2) from training at: /home/runner/work/EasyHybrid.jl/EasyHybrid.jl/docs/build
[ Info: Returning best model from epoch 199 of 300 epochs with best validation loss wrt mse: 0.014906584