From 8caab26e837067c928277947d3af17091c46f4a2 Mon Sep 17 00:00:00 2001 From: "Watal M. Iwasaki" Date: Wed, 2 Mar 2016 16:54:59 +0900 Subject: [PATCH] Add option for driver effects --- cell.cpp | 34 ++++++++++++++++++++++------------ cell.h | 19 +++++++++++-------- 2 files changed, 33 insertions(+), 20 deletions(-) diff --git a/cell.cpp b/cell.cpp index 335202a..cdf8e28 100644 --- a/cell.cpp +++ b/cell.cpp @@ -14,6 +14,7 @@ double Cell::MUTATION_RATE_ = 1e-1; double Cell::MUTATION_SIGMA_ = 0.0; +int Cell::DRIVER_EFFECTS_ = 1; double Cell::DRIVER_FRACTION_ = 0.0; double Cell::BIRTH_RATE_ = 1.0; double Cell::DEATH_RATE_ = 0.0; @@ -27,14 +28,16 @@ std::vector Cell::MUTANT_IDS_; //! Program options /*! @return Program options description - Command line option | Symbol | Variable - --------------------| ------------ | ------------------------- - `-u,--mutation` | \f$\mu\f$ | Cell::MUTATION_RATE_ - `-s,--sigma` | \f$\sigma\f$ | Cell::MUTATION_SIGMA_ - `-b,--birth` | \f$b\f$ | Cell::BIRTH_RATE_ - `-d,--death` | \f$d\f$ | Cell::DEATH_RATE_ - `-m,--migration` | \f$m\f$ | Cell::MIGRATION_RATE_ - `-k,--shape` | \f$k\f$ | Cell::GAMMA_SHAPE_ + Command line option | Symbol | Variable + --------------------| -------------- | ------------------------- + `-u,--mutation` | \f$\mu\f$ | Cell::MUTATION_RATE_ + `-s,--sigma` | \f$\sigma\f$ | Cell::MUTATION_SIGMA_ + `-e,--effect` | | Cell::DRIVER_EFFECTS_ + `-f,--fraction` | | Cell::DRIVER_FRACTION_ + `-b,--birth` | \f$\beta_0\f$ | Cell::BIRTH_RATE_ + `-d,--death` | \f$\delta_0\f$ | Cell::DEATH_RATE_ + `-m,--migration` | \f$m\f$ | Cell::MIGRATION_RATE_ + `-k,--shape` | \f$k\f$ | Cell::GAMMA_SHAPE_ */ boost::program_options::options_description& Cell::opt_description() { namespace po = boost::program_options; @@ -42,7 +45,8 @@ boost::program_options::options_description& Cell::opt_description() { desc.add_options() ("mutation,u", po::value(&MUTATION_RATE_)->default_value(MUTATION_RATE_)) ("sigma,s", po::value(&MUTATION_SIGMA_)->default_value(MUTATION_SIGMA_)) - ("driver,f", po::value(&DRIVER_FRACTION_)->default_value(DRIVER_FRACTION_)) + ("effect,e", po::value(&DRIVER_EFFECTS_)->default_value(DRIVER_EFFECTS_)) + ("fraction,f", po::value(&DRIVER_FRACTION_)->default_value(DRIVER_FRACTION_)) ("birth,b", po::value(&BIRTH_RATE_)->default_value(BIRTH_RATE_)) ("death,d", po::value(&DEATH_RATE_)->default_value(DEATH_RATE_)) ("migration,m", po::value(&MIGRATION_RATE_)->default_value(MIGRATION_RATE_)) @@ -65,7 +69,12 @@ void Cell::mutate() { sites_.push_back(MUTATION_EFFECTS_.size()); MUTATION_EFFECTS_.push_back(effect); MUTANT_IDS_.push_back(id()); - fitness_ *= (effect += 1.0); + if (DRIVER_EFFECTS_ & 0b01) { + birth_rate_ *= (effect += 1.0); + } + if (DRIVER_EFFECTS_ & 0b10) { + death_rate_ *= (1.0 - effect); + } } double Cell::delta_time() { @@ -110,7 +119,8 @@ std::string Cell::header(const size_t dimensions, const std::string& sep) { << "birth" << sep << "death" << sep; std::vector axes{"x", "y", "z"}; axes.resize(dimensions); - wtl::ost_join(oss, axes, sep) << sep << "sites" << sep << "fitness\n"; + wtl::ost_join(oss, axes, sep) << sep << "sites" << sep + << "beta" << sep << "delta\n"; return oss.str(); } @@ -119,7 +129,7 @@ std::ostream& Cell::write(std::ostream& ost, const std::string& sep) const { << time_of_birth_ << sep << time_of_death_ << sep; wtl::ost_join(ost, coord(), sep) << sep; wtl::ost_join(ost, sites(), ":") << sep - << fitness() << "\n"; + << birth_rate_ << death_rate_ << "\n"; return ost; } diff --git a/cell.h b/cell.h index 49f1af7..f6a5466 100644 --- a/cell.h +++ b/cell.h @@ -34,7 +34,8 @@ class Cell { Cell(const std::vector& v): coord_(v), id_(++ID_TAIL_), ancestor_(id_) {} //! Copy constructor Cell(const Cell& other): - coord_(other.coord_), sites_(other.sites_), fitness_(other.fitness_), + coord_(other.coord_), sites_(other.sites_), + birth_rate_(other.birth_rate_), death_rate_(other.death_rate_), id_(++ID_TAIL_), mother_(other.id_), ancestor_(other.ancestor_) {} //! Destructor ~Cell() {--ID_TAIL_;} @@ -71,8 +72,6 @@ class Cell { //! convert site positions to 01 vector std::vector haplotype(std::vector segsites) const; - double fitness() const {return fitness_;} - static const std::vector& MUTATION_EFFECTS() {return MUTATION_EFFECTS_;} static const std::vector& MUTANT_IDS() {return MUTANT_IDS_;} @@ -86,9 +85,9 @@ class Cell { private: //! finite per capita rates - double birth_rate() const {return BIRTH_RATE_ * fitness_;} - double death_rate() const {return DEATH_RATE_;} - double growth_rate() const {return 1.0 + birth_rate() - death_rate();} + double birth_rate() const {return birth_rate_;} + double death_rate() const {return death_rate_;} + double growth_rate() const {return 1.0 + birth_rate_ - death_rate_;} //! instantaneous rate for time increment double instantaneous_event_rate() const { const double lambda = growth_rate(); @@ -100,6 +99,9 @@ class Cell { static double MUTATION_SIGMA_; + //! 1: birth, 2: death, 3: both + static int DRIVER_EFFECTS_; + static double DRIVER_FRACTION_; static double BIRTH_RATE_; @@ -120,8 +122,9 @@ class Cell { std::vector coord_; //! Mutated sites (infinite-site model) std::vector sites_; - //! Sum of mutation effects + 1 - double fitness_ = 1.0; + + double birth_rate_ = BIRTH_RATE_; + double death_rate_ = DEATH_RATE_; //! Extra data size_t id_ = 0;