Commit dd371b96 authored by Jorn Bruggeman's avatar Jorn Bruggeman
Browse files

added model based on passive fluxes

parent 9f82e609
Strain Duplicate pCO2 Mol_Corg_per_cell Mol_Cinorg_per_cell Division_rate C_org_assim_rate C_inorg_assim_rate C_tot_assim_rate PICPOC Rs Rc d13Ccal d18Ocal d13Corg d13CDIC Sp_grad Cell_grad UAC Ll Th Volume CO2all CO2surf PH TEMP HCO3 CO3 CO2 dHCO3 dCO3 dCO2
1211 A 165.4136229 1.50E-12 1.86E-12 1.032684831 1.79E-17 2.22E-17 4.01E-17 1.242382271 3.738727584 2.929272973 -4.022 -5.512 -23.94892905 -3.8 0.18313595 0.110507352 0.5475 3.60742784 0.200105182 2.604075909 4.79E-06 4.00E-06 8.2 15 0.000871216 0.000102595 4.79E-06 -3.47694928 -6.07276992 -13.65383027
1211 A 305.7255676 9.65E-13 2.27E-12 0.995918826 1.11E-17 2.62E-17 3.73E-17 2.355319149 3.735819559 2.806797521 -3.364 -5.328 -23.73945704 -3.8 0.140128006 0.076926963 0.5475 3.570784442 0.195665394 2.49483191 8.85E-06 8.09E-06 8.2 15 0.001610224 0.000189621 8.85E-06 -3.47694928 -6.07276992 -13.65383027
1211 A 403.5577493 3.53E-12 2.99E-12 NA NA NA NA 0.84592145 3.582157738 2.835066964 -4.332 -5.385 -27.25118839 -3.8 0.111704379 0.095933228 0.5475 3.620839334 0.165822711 2.17401492 1.17E-05 NA 8.2 15 0.002125495 0.000250299 1.17E-05 -3.47694928 -6.07276992 -13.65383027
1211 A 1095.462981 1.84E-12 2.54E-12 0.988744303 2.10E-17 2.90E-17 5.00E-17 1.381974249 3.702702599 2.839082569 -3.558 -5.284 -26.66255956 -3.8 0.127999828 0.078048205 0.5475 3.481550107 0.155397237 1.883599616 3.17E-05 3.07E-05 8.2 15 0.005769686 0.00067944 3.17E-05 -3.47694928 -6.07276992 -13.65383027
1216 A 165.4136229 1.10E-12 8.77E-13 1.145219812 1.46E-17 1.16E-17 2.62E-17 0.795305164 2.63512502 2.221700311 -2.032 -3.909 -26.91602973 -3.8 0.083990951 0.05635936 0.585 2.507681002 0.084583616 0.531901024 4.79E-06 4.11E-06 8.2 15 0.000871216 0.000102595 4.79E-06 -3.47694928 -6.07276992 -13.65383027
1216 A 305.7255676 1.08E-12 9.13E-13 0.887152447 1.10E-17 9.38E-18 2.04E-17 0.849315068 2.619410519 2.200155055 -1.744 -4.724 -27.62242011 -3.8 0.070340187 0.045900334 0.585 2.420996593 0.087233961 0.511297831 8.85E-06 8.31E-06 8.2 15 0.001610224 0.000189621 8.85E-06 -3.47694928 -6.07276992 -13.65383027
1216 A 403.5577493 1.32E-12 1.03E-12 1.05398483 1.61E-17 1.25E-17 2.86E-17 0.779081134 2.570574348 2.125049331 -2.659 -4.806 -27.65654176 -3.8 0.067338702 0.05572572 0.585 2.476825507 0.079898802 0.490152352 1.17E-05 1.09E-05 8.2 15 0.002125495 0.000250299 1.17E-05 -3.47694928 -6.07276992 -13.65383027
1216 A 1095.462981 9.85E-13 8.07E-13 0.993215468 1.13E-17 9.28E-18 2.06E-17 0.819219791 2.51698967 NA -2.854 -4.901 -28.48871137 -3.8 0.054414812 0.133891108 0.585 2.404969811 0.084869003 0.490872109 3.17E-05 NA 8.2 15 0.005769686 0.00067944 3.17E-05 -3.47694928 -6.07276992 -13.65383027
1256 A 165.4136229 1.08E-12 4.56E-13 1.03619967 1.30E-17 5.47E-18 1.84E-17 0.422193878 2.463713768 2.200618961 -1.698 -3.114 -24.7110476 -3.8 0.061819418 0.053238509 0.39625 2.297223928 0.062846553 0.331656206 4.79E-06 4.31E-06 8.2 15 0.000871216 0.000102595 4.79E-06 -3.47694928 -6.07276992 -13.65383027
1256 A 305.7255676 9.29E-13 6.36E-13 1.037681473 1.12E-17 7.64E-18 1.88E-17 0.684254606 2.486266472 2.17626304 -2.06 -4.714 -23.8814697 -3.8 0.0620002 0.045869738 0.39625 2.297530537 0.061040876 0.32221321 8.85E-06 8.35E-06 8.2 15 0.001610224 0.000189621 8.85E-06 -3.47694928 -6.07276992 -13.65383027
1256 A 403.5577493 1.50E-12 7.29E-13 0.933446505 1.62E-17 7.88E-18 2.40E-17 0.487357775 2.518081819 2.184804412 -2.537 -4.265 -27.29414754 -3.8 0.058301141 0.044057445 0.39625 2.442814825 0.0620211 0.370101255 1.17E-05 1.10E-05 8.2 15 0.002125495 0.000250299 1.17E-05 -3.47694928 -6.07276992 -13.65383027
1256 A 1095.462981 1.27E-12 6.68E-13 1.05260628 1.54E-17 8.14E-18 2.36E-17 0.526610644 2.558617669 2.260123762 -2.63 -4.725 -25.54223372 -3.8 0.070943416 0.058468479 0.39625 2.441535478 0.061519428 0.366723185 3.17E-05 3.11E-05 8.2 15 0.005769686 0.00067944 3.17E-05 -3.47694928 -6.07276992 -13.65383027
1314 A 165.4136229 3.65E-12 3.67E-12 0.972907216 4.11E-17 4.13E-17 8.24E-17 1.005707763 4.317564672 3.563242948 -3.45 -5.182 -22.1507115 -3.8 0.309882932 0.204787307 0.5375 4.05390645 0.227735552 3.742641929 4.79E-06 3.47E-06 8.2 15 0.000871216 0.000102595 4.79E-06 -3.47694928 -6.07276992 -13.65383027
1314 A 305.7255676 2.80E-12 3.82E-12 0.959124138 3.11E-17 4.24E-17 7.35E-17 1.363636364 4.23412164 3.409303763 -3.111 -4.916 -23.25765479 -3.8 0.291706766 0.175173026 0.5375 4.03019479 0.197566488 3.20896777 8.85E-06 7.61E-06 8.2 15 0.001610224 0.000189621 8.85E-06 -3.47694928 -6.07276992 -13.65383027
1314 A 403.5577493 4.17E-12 4.13E-12 0.895550181 4.32E-17 4.28E-17 8.60E-17 0.990960452 4.233461669 3.471783978 -3.493 -5.407 -24.48508961 -3.8 0.269560467 0.177239536 0.5375 3.948556591 0.206192876 3.21477358 1.17E-05 1.03E-05 8.2 15 0.002125495 0.000250299 1.17E-05 -3.47694928 -6.07276992 -13.65383027
1314 A 1095.462981 3.47E-12 4.82E-12 0.774169283 3.11E-17 4.32E-17 7.43E-17 1.385964912 4.247272727 3.36 -3.418 -5.483 -23.83680575 -3.8 0.2027656 0.191647422 0.5375 3.869230818 0.240876285 3.60614613 3.17E-05 3.04E-05 8.2 15 0.005769686 0.00067944 3.17E-05 -3.47694928 -6.07276992 -13.65383027
1211 B 165.4136229 1.19E-12 2.18E-12 1.029514973 1.42E-17 2.59E-17 4.01E-17 1.831858407 3.723233696 2.897259964 -3.59 -4.651 -24.03115238 -3.8 0.179738896 0.106854393 0.5475 3.611192086 0.199672654 2.603872827 4.79E-06 3.99E-06 8.2 15 0.000871216 0.000102595 4.79E-06 -3.47694928 -6.07276992 -13.65383027
1211 B 305.7255676 1.26E-12 2.01E-12 0.984035036 1.43E-17 2.29E-17 3.72E-17 1.595744681 3.77107438 2.819820937 -3.312 -5.215 -23.30084147 -3.8 0.145705992 0.081362803 0.5475 3.601546815 0.222965682 2.892118959 8.85E-06 8.09E-06 8.2 15 0.001610224 0.000189621 8.85E-06 -3.47694928 -6.07276992 -13.65383027
1211 B 403.5577493 1.60E-12 2.80E-12 1.011051078 1.87E-17 3.28E-17 5.15E-17 1.752362949 3.650131867 2.82620915 -3.967 -5.461 -26.44811428 -3.8 0.141980826 0.100779015 0.5475 3.628122999 0.212335681 2.795033277 1.17E-05 1.06E-05 8.2 15 0.002125495 0.000250299 1.17E-05 -3.47694928 -6.07276992 -13.65383027
1211 B 1095.462981 1.73E-12 2.46E-12 0.870122223 1.74E-17 2.48E-17 4.22E-17 1.427098675 3.679064465 2.815974843 -3.658 -5.314 -26.09538938 -3.8 0.172029352 0.072714642 0.5475 3.655724651 0.212907659 2.845366663 3.17E-05 3.08E-05 8.2 15 0.005769686 0.00067944 3.17E-05 -3.47694928 -6.07276992 -13.65383027
1216 B 165.4136229 1.07E-12 8.62E-13 1.156833341 1.43E-17 1.15E-17 2.59E-17 0.806054872 2.626975633 2.220812231 -1.993 -4.004 NA -3.8 0.087687681 0.055960815 0.585 2.443181556 0.080460977 0.480282523 4.79E-06 4.12E-06 8.2 15 0.000871216 0.000102595 4.79E-06 -3.47694928 -6.07276992 -13.65383027
1216 B 305.7255676 8.25E-13 9.61E-13 1.006588126 9.61E-18 1.12E-17 2.08E-17 1.165876777 2.621626722 2.208011019 -1.783 -4.856 -27.568331 -3.8 0.078867155 0.046781934 0.585 2.420222721 0.08646041 0.506439954 8.85E-06 8.31E-06 8.2 15 0.001610224 0.000189621 8.85E-06 -3.47694928 -6.07276992 -13.65383027
1216 B 403.5577493 1.05E-12 9.19E-13 0.932911467 1.14E-17 9.93E-18 2.13E-17 0.872636816 2.547016591 2.109193816 -2.374 -4.704 -26.41416294 -3.8 0.068489283 0.055412371 0.585 2.533004491 0.08233226 0.528252982 1.17E-05 1.11E-05 8.2 15 0.002125495 0.000250299 1.17E-05 -3.47694928 -6.07276992 -13.65383027
1216 B 1095.462981 1.07E-12 8.23E-13 0.895004327 1.11E-17 8.53E-18 1.97E-17 0.766329347 2.523706851 2.161363636 -2.859 -4.888 -27.05656866 -3.8 0.08560033 0.044162966 0.585 2.39102619 0.074636541 0.42669757 3.17E-05 3.12E-05 8.2 15 0.005769686 0.00067944 3.17E-05 -3.47694928 -6.07276992 -13.65383027
1256 B 165.4136229 1.00E-12 5.22E-13 1.12681123 1.31E-17 6.81E-18 1.99E-17 0.52027972 2.469147382 2.186936653 -2.103 -4.275 -24.55481835 -3.8 0.060278422 0.050867811 0.39625 2.264036322 0.063579461 0.325899448 4.79E-06 4.26E-06 8.2 15 0.000871216 0.000102595 4.79E-06 -3.47694928 -6.07276992 -13.65383027
1256 B 305.7255676 9.30E-13 7.32E-13 1.064784496 1.15E-17 9.02E-18 2.05E-17 0.786713287 2.475812842 2.164993169 -2.203 -4.714 -24.2514019 -3.8 0.063170415 0.046029428 0.39625 2.414582893 0.063258001 0.368807464 8.85E-06 8.30E-06 8.2 15 0.001610224 0.000189621 8.85E-06 -3.47694928 -6.07276992 -13.65383027
1256 B 403.5577493 1.53E-12 5.92E-13 0.945030129 1.67E-17 6.47E-18 2.32E-17 0.388023952 2.525867294 2.197643678 -2.523 -4.342 -27.32915176 -3.8 0.063231836 0.046601604 0.39625 2.469310101 0.055492106 0.338362693 1.17E-05 1.11E-05 8.2 15 0.002125495 0.000250299 1.17E-05 -3.47694928 -6.07276992 -13.65383027
1256 B 1095.462981 1.09E-12 7.27E-13 1.106982496 1.40E-17 9.31E-18 2.33E-17 0.665869219 2.563658533 2.254795918 -2.83 -4.796 -26.06260933 -3.8 0.053232201 0.055955857 0.39625 2.406191448 0.063537251 0.367865262 3.17E-05 3.11E-05 8.2 15 0.005769686 0.00067944 3.17E-05 -3.47694928 -6.07276992 -13.65383027
1314 B 165.4136229 4.02E-12 3.79E-12 1.034414318 4.81E-17 4.54E-17 9.35E-17 0.942702703 4.340680697 3.579603762 -3.496 -5.202 -21.71748199 -3.8 0.329508537 0.220234373 0.5375 4.036990606 0.227328393 3.704837465 4.79E-06 3.29E-06 8.2 15 0.000871216 0.000102595 4.79E-06 -3.47694928 -6.07276992 -13.65383027
1314 B 305.7255676 3.69E-12 4.04E-12 0.833767886 3.56E-17 3.89E-17 7.45E-17 1.093902439 4.131315188 3.365994624 -3.088 -4.958 -23.56308534 -3.8 0.353899778 0.171627916 0.5375 4.146859677 0.184058918 3.165159101 8.85E-06 7.58E-06 8.2 15 0.001610224 0.000189621 8.85E-06 -3.47694928 -6.07276992 -13.65383027
1314 B 403.5577493 3.29E-12 3.44E-12 0.911535838 3.47E-17 3.63E-17 7.10E-17 1.044083527 4.205482821 3.444456893 -3.35 -5.436 -25.05344728 -3.8 0.256638569 0.167945951 0.5375 3.846225783 0.230846016 3.415009642 1.17E-05 1.05E-05 8.2 15 0.002125495 0.000250299 1.17E-05 -3.47694928 -6.07276992 -13.65383027
1314 B 1095.462981 2.93E-12 3.60E-12 1.008947336 3.42E-17 4.21E-17 7.63E-17 1.230566535 NA 3.340815603 -3.779 -5.385 -24.32641132 -3.8 0.252408591 0.197641664 0.5375 3.834027089 0.190314329 2.797575666 3.17E-05 3.04E-05 8.2 15 0.005769686 0.00067944 3.17E-05 -3.47694928 -6.07276992 -13.65383027
CBR A 261.95 2.10E-11 1.53E-11 0.970516646 2.36E-16 1.72E-16 4.07E-16 0.729837251 9.8945 8.230916902 -9.3285 -0.861 -21.8 -6.87 NA NA NA NA NA NA 9.10E-06 6.33E-06 8.13 18 0.001509318 0.000169621 9.10E-06 -6.573556156 -8.978714326 -16.37067084
CBR B 420.15 2.31E-11 9.48E-12 0.971759695 2.60E-16 1.07E-16 3.67E-16 0.410402275 10.3935 8.541285089 -9.1705 -0.9905 -21.7 -6.9 NA NA NA NA NA NA 1.46E-05 1.22E-05 8.13 18 0.002421543 0.000272139 1.46E-05 -6.603556156 -9.008641694 -16.40037498
CBR A 335.8 2.53E-11 1.96E-11 0.92798924 2.71E-16 2.10E-16 4.81E-16 0.774223175 9.9765 8.835200139 -9.02 -0.5835 -22.65 -7.035 NA NA NA NA NA NA 1.17E-05 8.66E-06 8.13 18 0.001940552 0.000218084 1.17E-05 -6.738556156 -9.143314849 -16.53404362
CBR B 364.1 2.48E-11 1.33E-11 0.928246921 2.67E-16 1.43E-16 4.09E-16 0.53419849 9.9355 8.779736978 -9.106 -0.823 -21.9 -6.995 NA NA NA NA NA NA 1.27E-05 1.00E-05 8.13 18 0.002098118 0.000235792 1.27E-05 -6.698556156 -9.103411692 -16.49443809
CBR A 497.25 2.16E-11 1.26E-11 0.967744744 2.41E-16 1.41E-16 3.83E-16 0.586022084 10.0375 8.318143286 -8.0715 -0.651 -23.3 -7.22 NA NA NA NA NA NA 1.73E-05 1.47E-05 8.13 18 0.002869363 0.000322466 1.73E-05 -6.923556156 -9.32786695 -16.71721916
CBR B 539.65 2.46E-11 1.51E-11 0.953425787 2.71E-16 1.66E-16 4.37E-16 0.612840556 9.827 8.743455409 -8.0355 -0.5645 -24.2 -6.96 NA NA NA NA NA NA 1.88E-05 1.60E-05 8.13 18 0.003109858 0.000349494 1.88E-05 -6.663556156 -9.068496429 -16.45978326
CBR A 1286.3 4.31E-11 2.85E-11 0.684578293 3.41E-16 2.26E-16 5.67E-16 0.661516633 9.13 10.82640502 -6.361 -0.06 -26.8 -7.05 NA NA NA NA NA NA 4.47E-05 4.18E-05 8.13 18 0.007413902 0.000833194 4.47E-05 -6.753556156 -9.158278533 -16.54889569
CBR B 1189.5 3.84E-11 2.65E-11 0.589261907 2.62E-16 1.81E-16 4.43E-16 0.69018455 9.025 10.36966142 -6.739 -0.934 -27.65 -7.06 NA NA NA NA NA NA 4.13E-05 3.89E-05 8.13 18 0.006849981 0.000769819 4.13E-05 -6.763556156 -9.168254322 -16.55879707
CBR B 1827.6 7.13E-11 3.56E-11 0.396089461 3.27E-16 1.63E-16 4.90E-16 0.499060974 NA 13.11893694 -6.2745 -0.2585 -27.95 -7.315 NA NA NA NA NA NA 6.35E-05 6.14E-05 8.13 18 0.01052376 0.001182688 6.35E-05 -7.018556156 -9.422636948 -16.81128227
import argparse
parser = argparse.ArgumentParser(description='Fit model to data (data+model defined in problem.py).')
parser.add_argument('--skipfit', action='store_true',help='Skip fitting and show result of current parameterization in model.py')
args = parser.parse_args()
import sys
sys.path.append('..')
import numpy
import problem
# Standard settings (minimum, mximum, log scale) for different parameters.
pardata = {'alpha':(1e-3,1e5,True),
'eps_FIX':(-30.,-10.,False),
'FB_vu_per_F_CAL':(0.,2.,False),
'k_f':(1e-1,1e1,True),
'c':(1e-3,1e1,True),
'n':(0.,5.,False),
'k_ei':(1e-5,1e-3,True),
'k_eiB':(1e-7,1e-3,True),
'k_plus_i_CA':(1e-6,3000.,True),
'Km':(0.01,0.2,False),
'f_chl':(0.0001,1.,True),
'f_cv':(0.0001,1.,True),
'T_a':(1e-5,1.,True),
'T_b':(1e-5,1.,True),
'k_plus_CA':(.1,1e8,True),
'pH_i':(7.,8.2,False),
'pH_x':(7.,8.2,False),
'pH_v':(7.,8.2,False),
}
# Create optimization problem (this will also create the model itself, as member "model")
pars = 'k_ei','T_a','T_b','eps_FIX','k_plus_CA' #,'pH_i','pH_x','pH_v'
prob = problem.Problem(pars)
prob.fit_FB_vu = False
# Create arrays with parameter bounds, initial estimates [taken from model.py], whether to log-transform.
# Optimization schemes used either initial estimates (SIMPLEX), or bounds (DIFFERENTIALEVOLUTION), but not both.
p_min = [pardata[par][0] for par in pars]
p_max = [pardata[par][1] for par in pars]
p_log = [pardata[par][2] for par in pars]
p_ini = [getattr(prob.model,par) for par in pars]
# If fitting FB_vu, add appropriate parameter bounds, initial estimates.
# Note we are fitting relative FB_vu between 0 and 1, which describes where
# the absolute FB_vu lies between the minimum and maximum possible values.
if prob.fit_FB_vu:
p_min += [0.]*len(problem.F_CALs)
p_max += [1.]*len(problem.F_CALs)
p_ini += [0.5]*len(problem.F_CALs)
p_log += [False]*len(problem.F_CALs)
if args.skipfit:
p = p_ini
else:
# Run the optimization
import optimize
method = optimize.SIMPLEX
method = (optimize.DIFFERENTIALEVOLUTION,optimize.SIMPLEX)
optimizer = optimize.Optimizer(prob)
p = optimizer.run(method,p_ini,par_min=p_min,par_max=p_max,logtransform=p_log)
print 'Optimization result:'
for par,value in zip(pars,p):
print ' %s = %s' % (par,value)
# Evaluate model at the optimum to ensure all parameters are set to their optimal value.
fitness,results = prob.evaluateFitness(p,return_results=True)
# Compute modelled isotopic compositions
deltaC_FIXs = numpy.empty_like(problem.d13Corgs)
deltaC_CALs = numpy.empty_like(problem.d13Ccals)
FB_ius = numpy.empty_like(problem.d13Ccals)
FB_xus = numpy.empty_like(problem.d13Ccals)
FB_vus = numpy.empty_like(problem.d13Ccals)
for i,res in enumerate(results):
deltaC_FIXs[i] = res['delta13C_org']
deltaC_CALs[i] = res['delta13C_cal']
FB_ius[i] = res['FB_iu']/(res['F_FIX']+res['F_CAL'])
FB_xus[i] = res['FB_xu']/res['F_FIX']
FB_vus[i] = res['FB_vu']/res['F_CAL']
# Plotting
import pylab
pylab.rcParams['mathtext.default']='regular'
def plotAll(x,y):
for istrain in problem.unique_strains:
valid = problem.strains==istrain
pylab.plot(x[valid],y[valid],'o',label=istrain)
pylab.figure(figsize=(15,8))
pylab.subplot(2,3,1)
plotAll(problem.d13Corgs-problem.d13CDICs,deltaC_FIXs-problem.d13CDICs)
pylab.xlabel(r'observed $\delta^{13}C_{org}$ - $\delta^{13}C_{DIC}$')
pylab.ylabel(r'modelled $\delta^{13}C_{org}$ - $\delta^{13}C_{DIC}$')
pylab.grid(True)
pylab.subplot(2,3,2)
plotAll(problem.pCO2s,deltaC_FIXs-problem.d13CDICs)
pylab.xlabel(r'$pCO_2$')
pylab.ylabel(r'modelled $\delta^{13}C_{org}$ - $\delta^{13}C_{DIC}$')
pylab.grid(True)
pylab.subplot(2,3,3)
plotAll(problem.pCO2s,problem.d13Corgs-problem.d13CDICs)
pylab.xlabel(r'$pCO_2$')
pylab.ylabel(r'observed $\delta^{13}C_{org}$ - $\delta^{13}C_{DIC}$')
#plotAll(problem.F_FIXs/problem.Cpercell/problem.pCO2s,deltaC_FIXs-problem.d13CDICs)
#pylab.xlabel(r'$F_{FIX}/V/pCO_2$')
#pylab.ylabel(r'modelled $\delta^{13}C_{org}$ - $\delta^{13}C_{DIC}$')
pylab.grid(True)
pylab.subplot(2,3,4)
plotAll(problem.d13Ccals-problem.d13CDICs,deltaC_CALs-problem.d13CDICs)
pylab.xlabel(r'observed $\delta^{13}C_{cal}$ - $\delta^{13}C_{DIC}$')
pylab.ylabel(r'modelled $\delta^{13}C_{cal}$ - $\delta^{13}C_{DIC}$')
pylab.grid(True)
pylab.subplot(2,3,5)
plotAll(problem.pCO2s,deltaC_CALs-problem.d13CDICs)
pylab.xlabel(r'$pCO_2$')
pylab.ylabel(r'modelled $\delta^{13}C_{cal}$ - $\delta^{13}C_{DIC}$')
pylab.grid(True)
pylab.legend(loc=1)
pylab.subplot(2,3,6)
plotAll(problem.pCO2s,problem.d13Ccals-problem.d13CDICs)
pylab.xlabel(r'$pCO_2$')
pylab.ylabel(r'observed $\delta^{13}C_{cal}$ - $\delta^{13}C_{DIC}$')
#plotAll(problem.F_FIXs/problem.Cpercell/problem.pCO2s,deltaC_CALs-problem.d13CDICs)
#pylab.xlabel(r'$F_{FIX}/V/pCO_2$')
#pylab.ylabel(r'modelled $\delta^{13}C_{cal}$ - $\delta^{13}C_{DIC}$')
pylab.grid(True)
pylab.savefig('match.png',dpi=150)
pylab.figure(figsize=(15,8))
pylab.subplot(3,1,1)
plotAll(problem.pCO2s,FB_ius)
pylab.xlabel(r'pCO2')
pylab.ylabel(r'FB_iu')
pylab.grid(True)
pylab.subplot(3,1,2)
plotAll(problem.pCO2s,FB_xus)
pylab.xlabel(r'pCO2')
pylab.ylabel(r'FB_xu')
pylab.grid(True)
pylab.subplot(3,1,3)
plotAll(problem.pCO2s,FB_vus)
pylab.xlabel(r'pCO2')
pylab.ylabel(r'FB_vu')
pylab.grid(True)
pylab.show()
# Stuff below not currently used, but aimed at computing actual energetic cost in W.
##
##print ssq,len(rs)
##saas
##FB_iu_bnds = numpy.linspace(0.,FB_iu*2,199)
##FB_xu_bnds = numpy.linspace(0.,FB_xu*2,200)
##
##FB_iu_ctr = 0.5*(FB_iu_bnds[:-1]+FB_iu_bnds[1:])
##FB_xu_ctr = 0.5*(FB_xu_bnds[:-1]+FB_xu_bnds[1:])
##
##FB_iu_ctr,FB_xu_ctr = numpy.meshgrid(FB_iu_ctr,FB_xu_ctr)
##
##cost = cost_function(FB_iu_ctr,FB_xu_ctr)
##print cost
##cost = numpy.ma.masked_invalid(cost)
##cost *= 8.31*(273.+15.) # convert from mol s-1 to W
##cost *= 1e12 # convert from W to pW
##mincost = cost.min()
##
##print 'Cost: %s - %s' % (cost.min(),cost.max())
##import pylab
###pylab.pcolormesh(FB_iu_bnds,FB_xu_bnds,cost,vmin=mincost,vmax=mincost*2)
##cf = pylab.contourf(FB_iu_ctr,FB_xu_ctr,cost,numpy.linspace(mincost,mincost+abs(mincost*2),20),extend='both')
##pylab.plot((FB_iu,),(FB_xu,),'o',mec='k',mfc='w')
##pylab.xlabel('FB_iu (mol s-1)')
##pylab.ylabel('FB_xu (mol s-1)')
##cb = pylab.colorbar(cf)
##cb.set_label('cost (pW)')
##pylab.savefig('FB_iu_FB_xu.png',dpi=300)
##pylab.show()
This diff is collapsed.
import scipy.optimize
import numpy
def oblateSpheroid(a,c):
e = numpy.sqrt(1-c*c/a/a)
S = 2.*numpy.pi*a*a + numpy.pi*c*c/e*numpy.log((1+e)/(1-e))
V = (4.*numpy.pi/3.)*a*a*c
return S,V
class Model(object):
def __init__(self):
# Cell membrane permeability (m/s)
self.k_ei = 0.0004 # m/s, Hopkinson et al. Table 1: 1.5 to 5.6 10^-2 cm/s, mean is 3.2 10^-2 cm/s
#self.k_plus_i_CA = 450 # s-1, Hopkinson et al. Table S3
#self.k_plus_x_CA = 2700 # s-1, Hopkinson et al. Table S3
#self.k_plus_v_CA = 0. # s-1, no CA
self.k_plus_CA = 4e7
self.pH_i = 7.
self.pH_x = 7.9
self.pH_v = 7.1
self.T_a = 0.0007
self.T_b = 0.012
self.f_x = 1.1
self.f_v = 0.6
self.S_x = 6.
self.S_v = 4.
#deltaC_e = deltaB_e - 9. # to ensure consistency with eps_cb-eps_bc
self.eps_CAL = 1.0 # per mille, Bolton & Stoll table S2: +1
self.eps_FIX = -15.0 # per mille, Bolton & Stoll table S2: -27
self.eps_cb_CA = -1.0 # per mille, Zeebe Table 3.2.5, p 184, 25 degrees Celsius!
self.eps_bc_CA = -10.0 # per mille, Zeebe Table 3.2.5, p 184, 25 degrees Celsius!
self.eps_cb_noCA = -13.0 # per mille, Zeebe Table 3.2.5, p 184, 24 degrees Celsius!
self.eps_bc_noCA = -22.0 # per mille, Zeebe Table 3.2.5, p 184, 24 degrees Celsius!
self.eps_cb4 = -11.
self.eps_bc4 = -20.
# To ensure difference is equal to equilibrium fractionation factor 10.3,
# which is used below to compute deltaC_e from deltaB_e.
#self.eps_cb *= 10.3/9.
#self.eps_bc *= 10.3/9.
def getBicarbonateFluxes(self,**info):
verbose = info.get('verbose',False)
if verbose:
for k,v in sorted(info.items(),cmp=lambda x,y: cmp(x[0],y[0])): print '%s = %s' % (k,v)
# Catch common unit errors
assert info['C_e']<1. and info['C_e']>0.001, 'Invalid C_e=%s. CO2 [aq] should lie between 0.001 and 1 mol m-3' % info['C_e']
assert info['r']<0.0001, 'Invalid r=%s. radius in m should be less than 0.0001' % info['r']
self.k_plus_i_CA = self.k_plus_CA
self.k_plus_x_CA = self.k_plus_CA
self.k_plus_v_CA = self.k_plus_CA
self.k_ex = self.k_ei
self.k_ev = self.k_ei
# Shortcuts to commonly used arguments
F_FIX = info['F_FIX']
F_CAL = info['F_CAL']
C_e = info['C_e']
B_e = info['B_e']
r = info['r']
H_i = 10.**-self.pH_i
H_v = 10.**-self.pH_v
H_x = 10.**-self.pH_x
SA = 4*numpy.pi*r**2 # cell surface area in m2, assuming spherical shape
V_cell = 4./3.*numpy.pi*r**3 # cell volume in m3, assuming spherical shape
if verbose: print 'division rate based on carbon estimated from volume = %.3f d-1' % (F_FIX/(20e3*V_cell)*86400,)
if verbose: print 'PIC:POC = %.3f' % (F_CAL/F_FIX,)
SA_x,V_x = oblateSpheroid(self.f_x*r,self.f_x*r/self.S_x)
SA_v,V_v = oblateSpheroid(self.f_v*r,self.f_v*r/self.S_v)
V_i = V_cell - V_v - V_x # cytosol volume (m3)
if verbose: print 'relative SA of chloroplast = %.3f, vesicle = %.3f' % (SA_x/SA,SA_v/SA)
if verbose: print 'relative volume of chloroplast = %.3f, vesicle = %.3f, cytosol = %.3f' % (V_x/V_cell,V_v/V_cell,V_i/V_cell)
Ut = (F_FIX+F_CAL)/(self.k_ei*C_e*SA + self.T_a*self.k_ei*B_e*SA)
P_Bi = self.k_ei*(self.T_a + self.T_b*Ut)
P_Bx = P_Bi
P_Bv = P_Bi
if verbose: print 'membrane permeability for HCO3 relative to CO2 = %.4f (utilization = %.3f)' % (P_Bi/self.k_ei,Ut)
# Compute k+ using temperature-dependent relation of Zeebe book, p 110, table 2.3.1
TK = info['temp'] + 273.15
k_plus_noCA = numpy.exp(1246.98 - 6.19e4/TK - 183.0*numpy.log(TK))
# Hydroxylation rate constant
R_igc = 8.3144621
k4_plus = 4.7e7 * numpy.exp(-23200/(R_igc*TK))
# Compute K1_star = [HCO3-]*[H+]/[CO2] = k_plus/k_minus
S = 35.
#pK1_star = 3670.7/TK - 62.008 + 9.7944*numpy.log(TK) - 0.0118*S + 0.000116*S*S # Millero p.664 (1995) using Mehrbach et al. data on seawater scale
pK1_star = 3633.86/TK - 61.2172 + 9.67770*numpy.log(TK) - 0.011555*S + 0.0001152*S*S # Lueker et al 2000 as per EPOCA recommendation
K1_star = 10.**(-pK1_star)
if verbose: print 'K1_star = %s' % K1_star
# Ionic product of water ([H+]*[OH-]) in M*M
# Source: EPOCA/OA recommendations Eq 1.10
K_W_star = numpy.exp(-13847.26/TK + 148.9652 - 23.6521*numpy.log(TK) + (118.67/TK - 5.977 + 1.0495*numpy.log(TK))*numpy.sqrt(S) - 0.01615*S)
OH_i = K_W_star/H_i
OH_x = K_W_star/H_x
OH_v = K_W_star/H_v
# Compute CO2-specific hydration rate (combine hydration rate constant with and without CA, and hydroxilation rate)
k_specHyd_i = k_plus_noCA + self.k_plus_i_CA + k4_plus*OH_i
k_specHyd_x = k_plus_noCA + self.k_plus_x_CA + k4_plus*OH_x
k_specHyd_v = k_plus_noCA + self.k_plus_v_CA + k4_plus*OH_v
if verbose: print 'cytosol CO2-specific hydration rate = %s s-1. CA contribution: %.2f %%' % (k_specHyd_i,100*self.k_plus_i_CA/k_specHyd_i)
if verbose: print 'chloroplast CO2-specific hydration rate = %s s-1. CA contribution: %.2f %%' % (k_specHyd_x,100*self.k_plus_x_CA/k_specHyd_x)
if verbose: print 'vesicle CO2-specific hydration rate = %s s-1. CA contribution: %.2f %%' % (k_specHyd_v,100*self.k_plus_v_CA/k_specHyd_v)
# Compute dehydration rate constants from hydration rate constants and equilibrium constant
# Units are s-1 M-1 - thus multiplication with [H+] in mol L-1 will give HCO3 dehydration rate in s-1
k_minus_noCA = k_plus_noCA/K1_star
k_minus_i_CA = self.k_plus_i_CA/K1_star
k_minus_x_CA = self.k_plus_x_CA/K1_star
k_minus_v_CA = self.k_plus_v_CA/K1_star
k4_minus = k4_plus*K_W_star/K1_star
k_specDeh_i = (k_minus_noCA + k_minus_i_CA)*H_i + k4_minus
k_specDeh_x = (k_minus_noCA + k_minus_x_CA)*H_x + k4_minus
k_specDeh_v = (k_minus_noCA + k_minus_v_CA)*H_v + k4_minus
if verbose: print 'cytosol HCO3-specific dehydration rate = %s s-1. CA contribution: %.1f %%' % (k_specDeh_i,k_minus_i_CA*H_i/k_specDeh_i*100)
if verbose: print 'chloroplast HCO3-specific dehydration rate = %s s-1. CA contribution: %.1f %%' % (k_specDeh_x,k_minus_x_CA*H_x/k_specDeh_x*100)
if verbose: print 'vesicle HCO3-specific dehydration rate = %s s-1. CA contribution: %.1f %%' % (k_specDeh_v,k_minus_v_CA*H_v/k_specDeh_v*100)
A = numpy.empty((6,6),dtype=float)
b = numpy.empty((6,),dtype=float)
A[0,:] = -self.k_ei*SA-self.k_ex*SA_x-self.k_ev*SA_v-k_specHyd_i*V_i,k_specDeh_i*V_i,self.k_ex*SA_x,0,self.k_ev*SA_v,0
A[1,:] = k_specHyd_i*V_i,-P_Bi*SA-P_Bx*SA_x-P_Bv*SA_v-k_specDeh_i*V_i,0,P_Bx*SA_x,0,P_Bv*SA_v
A[2,:] = self.k_ex*SA_x,0,-self.k_ex*SA_x-k_specHyd_x*V_x,k_specDeh_x*V_x,0,0
A[3,:] = 0,P_Bx*SA_x,k_specHyd_x*V_x,-P_Bx*SA_x-k_specDeh_x*V_x,0,0
A[4,:] = self.k_ev*SA_v,0,0,0,-self.k_ev*SA_v-k_specHyd_v*V_v,k_specDeh_v*V_v
A[5,:] = 0,P_Bv*SA_v,0,0,k_specHyd_v*V_v,-P_Bv*SA_v-k_specDeh_v*V_v
b[:] = -self.k_ei*SA*C_e,-P_Bi*SA*B_e,F_FIX,0,0,F_CAL
concs = numpy.linalg.solve(A,b)
if (concs<0).any(): return
C_i,B_i,C_x,B_x,C_v,B_v = concs
if verbose: print C_e,B_e,C_i,B_i,C_x,B_x,C_v,B_v
if verbose: print 'Actual B_i/C_i = %.1f. Equilibrium value = %.1f' % (B_i/C_i,K1_star/H_i)
if verbose: print 'Actual B_x/C_x = %.1f. Equilibrium value = %.1f' % (B_x/C_x,K1_star/H_x)
if verbose: print 'Actual B_v/C_v = %.1f. Equilibrium value = %.1f' % (B_v/C_v,K1_star/H_v)
A_iso = A*concs.reshape(1,-1)
A_iso[2,2] -= F_FIX
A_iso[5,5] -= F_CAL
# Passive CO2 fluxes from source concentration, membrane permeability, surface area.
FC_iu = self.k_ei*C_e*SA
FC_io = self.k_ei*C_i*SA
FC_xu = self.k_ex*C_i*SA_x
FC_xo = self.k_ex*C_x*SA_x
FC_vu = self.k_ev*C_i*SA_v
FC_vo = self.k_ev*C_v*SA_v
# Passive HCO3 fluxes from source concentration, membrane permeability, surface area.
FB_iu = P_Bi*B_e*SA
FB_io = P_Bi*B_i*SA
FB_xu = P_Bx*B_i*SA_x
FB_xo = P_Bx*B_x*SA_x
FB_vu = P_Bv*B_i*SA_v
FB_vo = P_Bv*B_v*SA_v
# CO2 hydration fluxes from CO2 concentration, rate constant, compartment volume.
FHyd_i = k_specHyd_i*C_i*V_i
FHyd_x = k_specHyd_x*C_x*V_x
FHyd_v = k_specHyd_v*C_v*V_v
# HCO3 dehydration fluxes from bicarbonate mass balance.
FDeh_i = k_specDeh_i*B_i*V_i
FDeh_x = k_specDeh_x*B_x*V_x
FDeh_v = k_specDeh_v*B_v*V_v
# Compute effective kinetic fractionation during hydration and dehydration
# from constants with and without CA.
eps_cb_i = (self.eps_cb_CA*self.k_plus_i_CA + self.eps_cb_noCA*k_plus_noCA + self.eps_cb4*k4_plus*OH_i)/k_specHyd_i
eps_cb_x = (self.eps_cb_CA*self.k_plus_x_CA + self.eps_cb_noCA*k_plus_noCA + self.eps_cb4*k4_plus*OH_x)/k_specHyd_x
eps_cb_v = (self.eps_cb_CA*self.k_plus_v_CA + self.eps_cb_noCA*k_plus_noCA + self.eps_cb4*k4_plus*OH_v)/k_specHyd_v
eps_bc_i = ((self.eps_bc_CA*k_minus_i_CA + self.eps_bc_noCA*k_minus_noCA)*H_i + self.eps_bc4*k4_minus)/k_specDeh_i
eps_bc_x = ((self.eps_bc_CA*k_minus_x_CA + self.eps_bc_noCA*k_minus_noCA)*H_x + self.eps_bc4*k4_minus)/k_specDeh_x
eps_bc_v = ((self.eps_bc_CA*k_minus_v_CA + self.eps_bc_noCA*k_minus_noCA)*H_v + self.eps_bc4*k4_minus)/k_specDeh_v
if verbose: print 'eps_cb_i = %s, eps_bc_i = %s' % (eps_cb_i,eps_bc_i)
if verbose: print 'eps_cb_x = %s, eps_bc_x = %s' % (eps_cb_x,eps_bc_x)
if verbose: print 'eps_cb_v = %s, eps_bc_v = %s' % (eps_cb_v,eps_bc_v)
# Find steady state isotopic composition by solving linear system.
# Matrix columns: deltaC_i,deltaB_i,deltaC_x,deltaB_x,deltaC_v,deltaB_v
# Matrix rows: time derivatives for C_i', B_i', C_x', B_x', C_v', B_v'
A[0,:] = (-FC_io-FC_xu-FC_vu-FHyd_i, FDeh_i, FC_xo, 0.0, FC_vo, 0.0)
A[1,:] = (FHyd_i, -FB_io-FB_xu-FB_vu-FDeh_i, 0.0, FB_xo, 0.0, FB_vo)
A[2,:] = (FC_xu, 0.0, -FC_xo-FHyd_x-F_FIX, FDeh_x, 0.0, 0.0)
A[3,:] = (0.0, FB_xu, FHyd_x, -FB_xo-FDeh_x, 0.0, 0.0)
A[4,:] = (FC_vu, 0.0, 0.0, 0.0, -FC_vo-FHyd_v, FDeh_v)
A[5,:] = (0.0, FB_vu, 0.0, 0.0, FHyd_v, -FB_vo-FDeh_v-F_CAL)
b[:] = (-eps_bc_i*FDeh_i + eps_cb_i*FHyd_i - info['deltaC_e']*FC_iu,
eps_bc_i*FDeh_i - eps_cb_i*FHyd_i - info['deltaB_e']*FB_iu,
-eps_bc_x*FDeh_x + eps_cb_x*FHyd_x + self.eps_FIX*F_FIX,
eps_bc_x*FDeh_x - eps_cb_x*FHyd_x,
-eps_bc_v*FDeh_v + eps_cb_v*FHyd_v,
eps_bc_v*FDeh_v - eps_cb_v*FHyd_v + self.eps_CAL*F_CAL)
deltaC_i,deltaB_i,deltaC_x,deltaB_x,deltaC_v,deltaB_v = numpy.linalg.solve(A,b)
if verbose: print deltaC_i,deltaB_i,deltaC_x,deltaB_x,deltaC_v,deltaB_v
# Compare alternative fomrualtons of matrix for isotope linear system
reldif = abs(A-A_iso)/abs(A).max()
assert (reldif<1e-13).all(),'Difference between isotope matrices too large: %s' % (reldif,)
# Compute delta13C of POC and PIC using isotopic composition of relevant pools and fractionation factors.
delta13C_org_mod = deltaC_x + self.eps_FIX
delta13C_cal_mod = deltaB_v + self.eps_CAL
if verbose: print delta13C_org_mod-info['deltaB_e']
result = info.copy()
result.update({'FB_iu':FB_iu,'FB_xu':FB_xu,'FB_vu':FB_vu,
'C_i':C_i,'C_x':C_x,'C_v':C_v,'C_e':C_e,
'deltaC_i':deltaC_i,'deltaB_i':deltaB_i,'deltaC_x':deltaC_x,
'deltaC_v':deltaC_v,'deltaB_v':deltaB_v,
'FC_io':FC_io,'FC_xo':FC_xo,'FC_vo':FC_vo,
'FC_iu':FC_iu,'FC_xu':FC_xu,'FC_vu':FC_vu,
'FHyd_i':FHyd_i,'FHyd_x':FHyd_x,'FHyd_v':FHyd_i,
'FDeh_i':FDeh_i,'FDeh_x':FDeh_x,'FDeh_v':FDeh_v,
'delta13C_org':delta13C_org_mod,'delta13C_cal':delta13C_cal_mod})
return result
import numpy
import model
m = model.Model()
# Model parameters (usually taken from best fit obtained earlier)
#m.k_ei = 0.004717235859
#m.T_a = 0.000814610301367
#m.T_b = 0.116871311168
#m.eps_FIX = -10.6106670595
#m.k_plus_CA = 2583.94260624
# Range of cell radius to try (micrometer)
rs = 1.,3.,9.
# Constant environment
temp = 15.
delta13C_HCO3 = 0. #0. #-3.5
pH = 8.2
S = 35.
eps_CO2_g = 0.0049*temp - 1.31 # From Zhang et al 1995
eps_HCO3_g = -0.114 *temp + 10.78 # From Zhang et al 1995
delta13C_CO2 = delta13C_HCO3 + eps_CO2_g - eps_HCO3_g
# Compute cytosolic K1_star = [HCO3-]*[H+]/[CO2] = k_plus_i/k_minus_i
TK = temp + 273.15
#pK1_star = 3670.7/TK - 62.008 + 9.7944*numpy.log(TK) - 0.0118*S + 0.000116*S*S # Millero p.664 (1995) using Mehrbach et al. data on seawater scale
pK1_star = 3633.86/TK - 61.2172 + 9.67770*numpy.log(TK) - 0.011555*S + 0.0001152*S*S # Lueker et al 2000 as per EPOCA recommendation
K1_star = 10.**(-pK1_star)
B_e_per_C_e = K1_star/10.**(-pH) # K1 = [HCO3][H]/[CO2] -> [HCO3] = K1*[CO2]/[H]
cmap = 'coolwarm_r'
cmap = 'jet'
import pylab
fig = pylab.figure()
for ir,r in enumerate(rs):
print 'Processing radius %s um' % r
# From radius in um to radius in m.
r = r*1e-6
# Volume of spherical cell (m3) from radius in m.
V = 4./3.*numpy.pi*r**3
# Estimate log10 of cellular organic carbon (mol) from cell volume (m3) using allometric relationship.
#log10C_org = 1.69904189333 + 0.836435899077*numpy.log10(V) # fit to Harry' data
#log10C_org = (numpy.log10(V)+18.)*0.991+numpy.log10(0.109) - numpy.log10(12e12) # Montagnes et al 1994: pg C from um3
#log10C_org = (numpy.log10(V)+18.)*0.875+numpy.log10(0.284) - numpy.log10(12e12) # Popp et al. 1998: pg C from um3
#C_org = 10.**log10C_org
# Cellular organic carbon in mol.
C_org = 20e3*V # Harry formulation 16-04-2015
# Use division rate of 1 per day and PIC:POC of 1.
div = .2
# PIC:POC range to evaluate
PIC_POCs = numpy.linspace(0.,2.5,100)
# pCO2 range to evaluate
pCO2s = numpy.linspace(100.,1200.,101)
#pCO2s = 367.,
delta13C_cals = numpy.ma.array(numpy.empty((pCO2s.size,PIC_POCs.size)),mask=True)
delta13C_orgs = numpy.ma.array(numpy.empty((pCO2s.size,PIC_POCs.size)),mask=True)
for i,pCO2 in enumerate(pCO2s):
for j,PIC_POC in enumerate(PIC_POCs):
# Compute fixation and calcification rates (mol/s)
F_FIX = C_org*div/86400.
F_CAL = F_FIX*PIC_POC
# Have model estimate bicarbonate fluxes and resulting delta13C values.
# Note pCO2 must be provided as fraction, i.e., ppm/1e6
C_e = 37.*pCO2*1e-6 # NB 37. from Zeebe book p 67 Fig 1.5.20, verify 37*367*1e-6 = 0.0136 mol m-3 = 13.6 uM
B_e = B_e_per_C_e*C_e
result = m.getBicarbonateFluxes(r=r,F_FIX=F_FIX,F_CAL=F_CAL,div=div,C_e=C_e,B_e=B_e,temp=temp,deltaB_e=delta13C_HCO3,deltaC_e=delta13C_CO2)
if result is not None:
delta13C_cals[i,j] = result['delta13C_cal']
delta13C_orgs[i,j] = result['delta13C_org']
pylab.figure(fig.number)
pylab.subplot(len(rs),2,ir*2+1)
cs = numpy.arange(numpy.floor(delta13C_cals.min()),numpy.ceil(delta13C_cals.max())+0.1,1.)
pc=pylab.contourf(PIC_POCs,pCO2s,delta13C_cals,cs,cmap=cmap)
pylab.contour(PIC_POCs,pCO2s,delta13C_cals,cs,colors='k',linewidths=.2,linestyles='-')
pylab.colorbar(pc)
if ir+1==len(rs): pylab.xlabel('PIC:POC (-)')
pylab.ylabel('pCO2 (ppm)')
pylab.subplot(len(rs),2,ir*2+2)
cs = numpy.arange(numpy.floor(delta13C_orgs.min()),numpy.ceil(delta13C_orgs.max())+0.1,1.)
pc=pylab.contourf(PIC_POCs,pCO2s,delta13C_orgs,cs,cmap=cmap)
pylab.contour(PIC_POCs,pCO2s,delta13C_orgs,cs,colors='k',linewidths=.2,linestyles='-')
pylab.colorbar(pc)
pylab.xlabel('PIC:POC (-)')
pylab.savefig('delta13C_trends.png',dpi=150)
pylab.show()
# This file combines the model with a set of observations
# to create an "optimization problem".
import sys,os
sys.path.append(r'C:\Users\Jorn\Documents')
import optimize
import numpy
import datafile
import model
df = datafile.DataFile()
# Filter dataset by requiring that the following columns do not have any missing data.
df.requireColumn('C_org_assim_rate')
df.requireColumn('C_inorg_assim_rate')
df.requireColumn('d13Ccal')
df.requireColumn('d13Corg')
df.requireColumn('Rc')
df.requireColumn('pCO2')
df.requireColumn('Strain')
df.requireColumn('Mol_Corg_per_cell')