WntCellCycleOdeSystem.cpp

00001 /*
00002 
00003 Copyright (C) University of Oxford, 2005-2010
00004 
00005 University of Oxford means the Chancellor, Masters and Scholars of the
00006 University of Oxford, having an administrative office at Wellington
00007 Square, Oxford OX1 2JD, UK.
00008 
00009 This file is part of Chaste.
00010 
00011 Chaste is free software: you can redistribute it and/or modify it
00012 under the terms of the GNU Lesser General Public License as published
00013 by the Free Software Foundation, either version 2.1 of the License, or
00014 (at your option) any later version.
00015 
00016 Chaste is distributed in the hope that it will be useful, but WITHOUT
00017 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
00018 FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
00019 License for more details. The offer of Chaste under the terms of the
00020 License is subject to the License being interpreted in accordance with
00021 English Law and subject to any action against the University of Oxford
00022 being under the jurisdiction of the English Courts.
00023 
00024 You should have received a copy of the GNU Lesser General Public License
00025 along with Chaste. If not, see <http://www.gnu.org/licenses/>.
00026 
00027 */
00028 #include "WntCellCycleOdeSystem.hpp"
00029 #include "CellwiseOdeSystemInformation.hpp"
00030 
00031 
00032 WntCellCycleOdeSystem::WntCellCycleOdeSystem(double wntLevel,
00033                                              boost::shared_ptr<AbstractCellMutationState> pMutationState,
00034                                              std::vector<double> stateVariables)
00035     : AbstractOdeSystem(9),
00036       mpMutationState(pMutationState),
00037       mWntLevel(wntLevel)
00038 {
00039     mpSystemInfo.reset(new CellwiseOdeSystemInformation<WntCellCycleOdeSystem>);
00040 
00055     Init(); // set up parameter values
00056 
00057     // Set up a Wnt signalling pathway in a steady state
00058     double destruction_level = ma5d/(ma4d*wntLevel+ma5d);
00059     double beta_cat_level_1 = -1.0;
00060     double beta_cat_level_2 = -1.0;
00061 
00062     if (!mpMutationState)
00063     {
00064         // No mutations specified
00065     }
00066     else if (mpMutationState && mpMutationState->IsType<ApcOneHitCellMutationState>())
00067     {
00068         // APC +/- : only half are active
00069         beta_cat_level_1 = 0.5*ma2d/(ma2d+0.5*ma3d*destruction_level);
00070         beta_cat_level_2 = 0.5*ma2d/(ma2d+0.5*ma3d*destruction_level);
00071     }
00072     else if (mpMutationState && mpMutationState->IsType<ApcTwoHitCellMutationState>())
00073     {
00074         // APC -/-
00075         destruction_level = 0.0; // no active destruction complex
00076         beta_cat_level_1 = 0.5; // fully active beta-catenin
00077         beta_cat_level_2 = 0.5; // fully active beta-catenin
00078     }
00079     else if (mpMutationState && mpMutationState->IsType<BetaCateninOneHitCellMutationState>())
00080     {
00081         // Beta-cat delta 45
00082         beta_cat_level_1 = 0.5*ma2d/(ma2d+ma3d*destruction_level);
00083         beta_cat_level_2 = 0.5;
00084     }
00085     else
00086     {
00087         // healthy cells
00088         beta_cat_level_1 = 0.5*ma2d/(ma2d+ma3d*destruction_level);
00089         beta_cat_level_2 = 0.5*ma2d/(ma2d+ma3d*destruction_level);
00090     }
00091 
00092     // Cell-specific initial conditions
00093     SetDefaultInitialCondition(5, destruction_level);
00094     SetDefaultInitialCondition(6, beta_cat_level_1);
00095     SetDefaultInitialCondition(7, beta_cat_level_2);
00096     SetDefaultInitialCondition(8, wntLevel);
00097 
00098     if (stateVariables != std::vector<double>())
00099     {
00100         SetStateVariables(stateVariables);
00101     }
00102 }
00103 
00104 void WntCellCycleOdeSystem::SetMutationState(boost::shared_ptr<AbstractCellMutationState> pMutationState)
00105 {
00106     mpMutationState = pMutationState;
00107 }
00108 
00109 WntCellCycleOdeSystem::~WntCellCycleOdeSystem()
00110 {
00111     // Do nothing
00112 }
00113 
00114 void WntCellCycleOdeSystem::Init()
00115 {
00116     // Initialise model parameter values
00117     // Swat (2004) Parameters
00118     double k1 = 1.0;
00119     double k2 = 1.6;
00120     double k3 = 0.05;
00121     double k16 = 0.4;
00122     double k34 = 0.04;
00123     double k43 = 0.01;
00124     double k61 = 0.3;
00125     double k23 = 0.3;
00126     double a = 0.04;
00127     double J11 = 0.5;
00128     double J12 = 5.0;
00129     double J61 = 5.0;
00130     double J62 = 8.0;
00131     double J13 = 0.002;
00132     double J63 = 2.0;
00133     double Km1 = 0.5;
00134     double Km2 = 4.0;
00135     double Km4 = 0.3;
00136     double kp = 0.05;
00137     double phi_pRb = 0.005;
00138     double phi_E2F1 = 0.1;
00139     double phi_CycDi = 0.023;
00140     double phi_CycDa = 0.03;
00141     double phi_pRbp = 0.06;
00142 
00143     // Mirams et al. parameter values
00144     double a1 = 0.423;
00145     double a2 = 2.57e-4;
00146     double a3 = 1.72;
00147     double a4 = 10.0;
00148     double a5 = 0.5;
00149     double WntMax = 10.0;
00150     double mitogenic_factorF = 6.0e-4;
00151     double APC_Total = 0.02;
00152 
00153     // Non-dimensionalise...
00154     mk2d = k2/(Km2*phi_E2F1);
00155     mk3d = k3*a1*mitogenic_factorF/(Km4*phi_E2F1*a2);
00156     mk34d = k34/phi_E2F1;
00157     mk43d = k43/phi_E2F1;
00158     mk23d = k23*Km2/(Km4*phi_E2F1);
00159     mad = a/Km2;
00160     mJ11d = J11*phi_E2F1/k1;
00161     mJ12d = J12*phi_E2F1/k1;
00162     mJ13d = J13*phi_E2F1/k1;
00163     mJ61d = J61*phi_E2F1/k1;
00164     mJ62d = J62*phi_E2F1/k1;
00165     mJ63d = J63*phi_E2F1/k1;
00166     mKm1d = Km1/Km2;
00167     mkpd = kp/(Km2*phi_E2F1);
00168     mphi_r = phi_pRb/phi_E2F1;
00169     mphi_i = phi_CycDi/phi_E2F1;
00170     mphi_j = phi_CycDa/phi_E2F1;
00171     mphi_p = phi_pRbp/phi_E2F1;
00172     ma2d = a2/phi_E2F1;
00173     ma3d = a3*APC_Total/phi_E2F1;
00174     ma4d = a4*WntMax/phi_E2F1;
00175     ma5d = a5/phi_E2F1;
00176     mk16d = k16*Km4/phi_E2F1;
00177     mk61d = k61/phi_E2F1;
00178     mPhiE2F1 = phi_E2F1;
00179 }
00180 
00181 void WntCellCycleOdeSystem::EvaluateYDerivatives(double time, const std::vector<double>& rY, std::vector<double>& rDY)
00182 {
00183     double r = rY[0];
00184     double e = rY[1];
00185     double i = rY[2];
00186     double j = rY[3];
00187     double p = rY[4];
00188     double c = rY[5];
00189     double b1 = rY[6];
00190     double b2 = rY[7];
00191     double wnt_level = rY[8];
00192 
00193     double dx1 = 0.0;
00194     double dx2 = 0.0;
00195     double dx3 = 0.0;
00196     double dx4 = 0.0;
00197     double dx5 = 0.0;
00198     double dx6 = 0.0;
00199     double dx7 = 0.0;
00200     double dx8 = 0.0;
00201 
00202     /*
00203      * The variables are
00204      * 1. r = pRb
00205      * 2. e = E2F1
00206      * 3. i = CycD (inactive)
00207      * 4. j = CycD (active)
00208      * 5. p = pRb-p
00209      * 6. c = APC (Active)
00210      * 7. b = Beta-Catenin
00211     */
00212 
00213     // Bit back-to-front, but work out the Wnt section first...
00214 
00215     // Mutations take effect by altering the level of beta-catenin
00216     if (!mpMutationState)
00217     {
00218         // No mutations specified
00219     }
00220     else if (mpMutationState->IsType<ApcOneHitCellMutationState>()) // APC +/-
00221     {
00222         dx6 = ma5d*(1.0-c) - ma4d*wnt_level*c;
00223         dx7 = ma2d*(0.5-b1) - 0.5*ma3d*b1*c;
00224         dx8 = ma2d*(0.5-b2) - 0.5*ma3d*b2*c;
00225     }
00226     else if (mpMutationState->IsType<ApcTwoHitCellMutationState>()) // APC -/-
00227     {
00228         dx6 = 0.0;
00229         dx7 = ma2d*(0.5-b1);
00230         dx8 = ma2d*(0.5-b2);
00231     }
00232     else if (mpMutationState->IsType<BetaCateninOneHitCellMutationState>()) // Beta-Cat D45
00233     {
00234         dx6 = ma5d*(1.0-c) - ma4d*wnt_level*c;
00235         dx7 = ma2d*(0.5-b1) - ma3d*b1*c;
00236         dx8 = ma2d*(0.5-b2);
00237     }
00238     else
00239     {
00240         // da
00241         dx6 = ma5d*(1.0-c) - ma4d*wnt_level*c;
00242         // db
00243         dx7 = ma2d*(0.5-b1) - ma3d*b1*c;
00244         dx8 = ma2d*(0.5-b2) - ma3d*b2*c;
00245     }
00246 
00247     // Now the cell cycle stuff...
00248 
00249     // dr
00250     dx1 = e/(mKm1d+e)*mJ11d/(mJ11d+r)*mJ61d/(mJ61d+p) - mk16d*r*j+mk61d*p-mphi_r*r;
00251     // de
00252     dx2 = mkpd+mk2d*(mad*mad+e*e)/(1+e*e)*mJ12d/(mJ12d+r)*mJ62d/(mJ62d+p) - e;
00253     // di
00254     dx3 = mk3d*(b1+b2) + mk23d*e*mJ13d/(mJ13d+r)*mJ63d/(mJ63d+p) + mk43d*j - mk34d*i*j/(1+j) - mphi_i*i;
00255     // dj
00256     dx4 = mk34d*i*j/(1+j) - (mk43d+mphi_j)*j;
00257     // dp
00258     dx5 = mk16d*r*j - mk61d*p - mphi_p*p;
00259 
00260     double factor = mPhiE2F1*60.0;  // convert non-dimensional d/dt s to d/dt in hours
00261 
00262     rDY[0] = dx1*factor;
00263     rDY[1] = dx2*factor;
00264     rDY[2] = dx3*factor;
00265     rDY[3] = dx4*factor;
00266     rDY[4] = dx5*factor;
00267     rDY[5] = dx6*factor;
00268     rDY[6] = dx7*factor; // beta-cat allele 1
00269     rDY[7] = dx8*factor; // beta-cat allele 2
00270     rDY[8] = 0.0; // do not change the Wnt level
00271 }
00272 
00273 const boost::shared_ptr<AbstractCellMutationState> WntCellCycleOdeSystem::GetMutationState() const
00274 {
00275     return mpMutationState;
00276 }
00277 
00278 bool WntCellCycleOdeSystem::CalculateStoppingEvent(double time, const std::vector<double>& rY)
00279 {
00280     double r = rY[0];
00281     double e = rY[1];
00282     double p = rY[4];
00283     double dY1 = mkpd+mk2d*(mad*mad+e*e)/(1+e*e)*mJ12d/(mJ12d+r)*mJ62d/(mJ62d+p) - e;
00284     double factor = mPhiE2F1*60.0;  // Convert non-dimensional d/dt s to d/dt in hours.
00285     dY1 = dY1*factor;
00286 
00287     assert(!std::isnan(rY[1]));
00288     assert(!std::isnan(dY1));
00289     return (rY[1] > 1.0 && dY1 > 0.0);
00290 }
00291 
00292 double WntCellCycleOdeSystem::CalculateRootFunction(double time, const std::vector<double>& rY)
00293 {
00294     return rY[1] - 1.0;
00295 }
00296 
00297 template<>
00298 void CellwiseOdeSystemInformation<WntCellCycleOdeSystem>::Initialise()
00299 {
00300     this->mVariableNames.push_back("pRb");
00301     this->mVariableUnits.push_back("non_dim");
00302     this->mInitialConditions.push_back(7.357000000000000e-01);
00303 
00304     this->mVariableNames.push_back("E2F1");
00305     this->mVariableUnits.push_back("non_dim");
00306     this->mInitialConditions.push_back(1.713000000000000e-01);
00307 
00308     this->mVariableNames.push_back("CycD_i");
00309     this->mVariableUnits.push_back("non_dim");
00310     this->mInitialConditions.push_back(6.900000000000001e-02);
00311 
00312     this->mVariableNames.push_back("CycD_a");
00313     this->mVariableUnits.push_back("non_dim");
00314     this->mInitialConditions.push_back(3.333333333333334e-03);
00315 
00316     this->mVariableNames.push_back("pRb_p");
00317     this->mVariableUnits.push_back("non_dim");
00318     this->mInitialConditions.push_back(1.000000000000000e-04);
00319 
00320     this->mVariableNames.push_back("APC");
00321     this->mVariableUnits.push_back("non_dim");
00322     this->mInitialConditions.push_back(NAN); // will be filled in later
00323 
00324     this->mVariableNames.push_back("Beta_Cat1");
00325     this->mVariableUnits.push_back("non_dim");
00326     this->mInitialConditions.push_back(NAN); // will be filled in later
00327 
00328     this->mVariableNames.push_back("Beta_Cat2");
00329     this->mVariableUnits.push_back("non_dim");
00330     this->mInitialConditions.push_back(NAN); // will be filled in later
00331 
00332     this->mVariableNames.push_back("Wnt");
00333     this->mVariableUnits.push_back("non_dim");
00334     this->mInitialConditions.push_back(NAN); // will be filled in later
00335 
00336     this->mInitialised = true;
00337 }
00338 
00339 double WntCellCycleOdeSystem::GetWntLevel() const
00340 {
00341     return mWntLevel;
00342 }
00343 
00344 // Serialization for Boost >= 1.36
00345 #include "SerializationExportWrapperForCpp.hpp"
00346 CHASTE_CLASS_EXPORT(WntCellCycleOdeSystem)

Generated on Mon Nov 1 12:35:16 2010 for Chaste by  doxygen 1.5.5