Menyanthes File

In [1]:

    

# First perform the necessary imports
import matplotlib.pyplot as plt
import pastas as ps

%matplotlib notebook


    

In [2]:

    

# how to use it?
fname = '../data/MenyanthesTest.men'
meny = ps.read.MenyData(fname)

# plot some series
f1, axarr = plt.subplots(len(meny.IN)+1, sharex=True)
oseries = meny.H['Obsevation well']["values"]
oseries.plot(ax=axarr[0])
axarr[0].set_title(meny.H['Obsevation well']["Name"])
for i, val in enumerate(meny.IN.items()):
    name, data = val
    data["values"].plot(ax=axarr[i+1])
    axarr[i+1].set_title(name)
plt.tight_layout(pad=0)
plt.show()


    

In [3]:

    

# Create the time series model
ml = ps.Model(oseries)

# Add precipitation
IN = meny.IN['Precipitation']['values']
IN.index = IN.index.round("D")
IN2 = meny.IN['Evaporation']['values']
IN2.index = IN2.index.round("D")
ts = ps.StressModel2([IN, IN2], ps.Gamma, 'Recharge')
ml.add_stressmodel(ts)

# Add well extraction 1
# IN = meny.IN['Extraction 1']
# # extraction amount counts for the previous month
# ts = ps.StressModel(IN['values'], ps.Hantush, 'Extraction_1', up=False,
#                     settings="well")
# ml.add_stressmodel(ts)

# Add well extraction 2
IN = meny.IN['Extraction 2']
# extraction amount counts for the previous month
ts = ps.StressModel(IN['values'], ps.Hantush, 'Extraction_2', up=False,
                    settings="well")
ml.add_stressmodel(ts)

# Add well extraction 3
IN = meny.IN['Extraction 3']
# extraction amount counts for the previous month
ts = ps.StressModel(IN['values'], ps.Hantush, 'Extraction_3', up=False,
                    settings="well")
ml.add_stressmodel(ts)

# Solve the model (can take around 20 seconds..)
ml.solve()


    

    




INFO: Cannot determine frequency of series None
INFO: Inferred frequency from time series None: freq=D 
INFO: Inferred frequency from time series None: freq=D 
INFO: Cannot determine frequency of series None
INFO: Time Series None: values of stress were transformedto daily values (frequency not altered) with: divide
INFO: Time Series None: values of stress were transformedto daily values (frequency not altered) with: divide
INFO: Time Series None was sampled down to freq D with method timestep_weighted_resample
INFO: Cannot determine frequency of series None
INFO: Time Series None: values of stress were transformedto daily values (frequency not altered) with: divide
INFO: Time Series None: values of stress were transformedto daily values (frequency not altered) with: divide
INFO: Time Series None was sampled down to freq D with method timestep_weighted_resample
INFO: Time Series None: values of stress were transformedto daily values (frequency not altered) with: divide
INFO: Time Series None was sampled down to freq D with method timestep_weighted_resample
INFO: Time Series None: values of stress were transformedto daily values (frequency not altered) with: divide
INFO: Time Series None was sampled down to freq D with method timestep_weighted_resample
INFO: There are observations between the simulation timesteps. Linear interpolation between simulated values is used.
/Users/Raoul/Projects/pastas/pastas/pastas/solver.py:118: RuntimeWarning: invalid value encountered in double_scalars
  pcor[i, j] = pcov[i, j] / np.sqrt(pcov[i, i] * pcov[j, j])






    




Model Results Observations                Fit Statistics
========================================================
nfev     8                      EVP                93.81
nobs     2843                   R2                  0.94
noise    True                   RMSE                0.22
tmin     1960-04-28 00:00:00    AIC                14.09
tmax     2015-06-29 00:00:00    BIC                85.53
freq     D                      ___                     
warmup   3650 days 00:00:00     ___                     
solver   LeastSquares           ___                     

Parameters (12 were optimized)
========================================================
                      optimal   stderr     initial  vary
Recharge_A        1769.558531   ±3.28%  210.498526  True
Recharge_n           1.041631   ±1.52%    1.000000  True
Recharge_a         718.996058   ±5.93%   10.000000  True
Recharge_f          -1.048179   ±3.86%   -1.000000  True
Extraction_2_A      -0.003745   ±0.95%   -0.002620  True
Extraction_2_rho     0.182700   ±9.90%    1.000000  True
Extraction_2_cS   1826.220540  ±10.14%  100.000000  True
Extraction_3_A      -0.001004   ±3.76%   -0.005146  True
Extraction_3_rho     4.050843  ±42.82%    1.000000  True
Extraction_3_cS     45.029454  ±43.29%  100.000000  True
constant_d          10.184589   ±1.12%    8.557530  True
noise_alpha          0.013454   ±0.00%    1.000000  True

Parameter correlations |rho| > 0.5
========================================================
Recharge_A       Recharge_a        0.73
                 Extraction_3_A    0.51
Recharge_n       Recharge_a       -0.69
                 Extraction_2_rho -0.54
Recharge_a       Extraction_2_cS  -0.55
Recharge_f       constant_d       -0.95
Extraction_2_A   Extraction_3_A   -0.70
Extraction_2_rho Extraction_2_cS  -0.87
Extraction_3_rho Extraction_3_cS  -0.97

In [4]:

    

ax = ml.plots.decomposition(ytick_base=1.)
ax[0].set_title('Observations vs simulation')
ax[0].legend()
ax[0].figure.tight_layout(pad=0)


    

In [ ]: