Skip to content

Commit

Permalink
Fix handling of "super lethal" mutations. (#1319)
Browse files Browse the repository at this point in the history
* for multiplicative fitness, return 0 as soon as w <= 0.
* extract C++ back end to locations that enable C++
  testing in cpptests/
  • Loading branch information
molpopgen authored Jul 7, 2024
1 parent bee37c2 commit 0d36339
Show file tree
Hide file tree
Showing 10 changed files with 687 additions and 312 deletions.
82 changes: 5 additions & 77 deletions cpp/genetic_values/Additive.cc
Original file line number Diff line number Diff line change
Expand Up @@ -21,97 +21,25 @@
#include <stdexcept>

#include <fwdpy11/types/Mutation.hpp>
#include <fwdpy11/genetic_values/fwdpp_wrappers/fwdpp_genetic_value.hpp>
#include <fwdpy11/genetic_values/DiploidAdditive.hpp>
#include <fwdpy11/genetic_value_to_fitness/GeneticValueIsTrait.hpp>
#include <pybind11/pybind11.h>

namespace
{
struct single_deme_additive_het
{
inline void
operator()(double& d, const fwdpy11::Mutation& m) const
{
d += m.s * m.h;
}
};

struct multi_deme_additive_het
{
inline void
operator()(const std::size_t deme, double& d, const fwdpy11::Mutation& m) const
{
d += m.esizes[deme] * m.heffects[deme];
}
};

struct single_deme_additive_hom
{
double scaling;
single_deme_additive_hom(double s) : scaling(s)
{
}

inline void
operator()(double& d, const fwdpy11::Mutation& m) const
{
d += scaling * m.s;
}
};

struct multi_deme_additive_hom
{
double scaling;
multi_deme_additive_hom(double s) : scaling(s)
{
}

inline void
operator()(const std::size_t deme, double& d, const fwdpy11::Mutation& m) const
{
d += scaling * m.esizes[deme];
}
};

struct final_additive_trait
{
inline double
operator()(double d) const
{
return d;
}
};

struct final_additive_fitness
{
inline double
operator()(double d) const
{
return std::max(0.0, 1. + d);
}
};

using DiploidAdditive = fwdpy11::stateless_site_dependent_genetic_value_wrapper<
single_deme_additive_het, single_deme_additive_hom, multi_deme_additive_het,
multi_deme_additive_hom, 0>;
}

namespace py = pybind11;

void
init_Additive(py::module& m)
{
py::class_<DiploidAdditive, fwdpy11::DiploidGeneticValue>(m, "_ll_Additive")
py::class_<fwdpy11::DiploidAdditive, fwdpy11::DiploidGeneticValue>(m, "_ll_Additive")
.def(py::init([](double scaling,
const fwdpy11::GeneticValueIsTrait* gvalue_to_fitness,
const fwdpy11::GeneticValueNoise* noise, std::size_t ndemes) {
if (gvalue_to_fitness != nullptr)
{
return DiploidAdditive(ndemes, scaling, final_additive_trait(),
gvalue_to_fitness, noise);
return fwdpy11::additive_trait_model(ndemes, scaling,
gvalue_to_fitness, noise);
}
return DiploidAdditive(ndemes, scaling, final_additive_fitness(),
gvalue_to_fitness, noise);
return fwdpy11::additive_fitness_model(ndemes, scaling, noise);
}),
py::arg("scaling"), py::arg("gvalue_to_fitness"), py::arg("noise"),
py::arg("ndemes"));
Expand Down
88 changes: 6 additions & 82 deletions cpp/genetic_values/Multiplicative.cc
Original file line number Diff line number Diff line change
Expand Up @@ -17,104 +17,28 @@
// along with fwdpy11. If not, see <http://www.gnu.org/licenses/>.
//
#include <functional>
#include <stdexcept>

#include <fwdpy11/genetic_value_to_fitness/GeneticValueIsTrait.hpp>
#include <fwdpy11/genetic_values/fwdpp_wrappers/fwdpp_genetic_value.hpp>
#include <fwdpy11/genetic_values/DiploidMultiplicative.hpp>
#include <fwdpy11/types/Mutation.hpp>
#include <pybind11/pybind11.h>

namespace py = pybind11;

namespace
{
struct single_deme_multiplicative_het
{
inline void
operator()(double& d, const fwdpy11::Mutation& m) const
{
d *= (1. + m.s * m.h);
}
};

struct multi_deme_multiplicative_het
{
inline void
operator()(const std::size_t deme, double& d, const fwdpy11::Mutation& m) const
{
d *= (1. + m.esizes[deme] * m.heffects[deme]);
}
};

struct single_deme_multiplicative_hom
{
double scaling;
single_deme_multiplicative_hom(double s) : scaling(s)
{
}

inline void
operator()(double& d, const fwdpy11::Mutation& m) const
{
d *= (1. + scaling * m.s);
}
};

struct multi_deme_multiplicative_hom
{
double scaling;
multi_deme_multiplicative_hom(double s) : scaling(s)
{
}

inline void
operator()(const std::size_t deme, double& d, const fwdpy11::Mutation& m) const
{
d *= (1. + scaling * m.esizes[deme]);
}
};

struct final_multiplicative_trait
{
inline double
operator()(double d) const
{
return d - 1.0;
}
};

struct final_multiplicative_fitness
{
inline double
operator()(double d) const
{
return std::max(0.0, d);
}
};

using DiploidMultiplicative
= fwdpy11::stateless_site_dependent_genetic_value_wrapper<
single_deme_multiplicative_het, single_deme_multiplicative_hom,
multi_deme_multiplicative_het, multi_deme_multiplicative_hom, 1>;
}

void
init_Multiplicative(py::module& m)
{
py::class_<DiploidMultiplicative, fwdpy11::DiploidGeneticValue>(m,
"_ll_Multiplicative")
py::class_<fwdpy11::DiploidMultiplicative, fwdpy11::DiploidGeneticValue>(
m, "_ll_Multiplicative")
.def(py::init([](double scaling,
const fwdpy11::GeneticValueIsTrait* gvalue_to_fitness,
const fwdpy11::GeneticValueNoise* noise, std::size_t ndemes) {
if (gvalue_to_fitness != nullptr)
{
return DiploidMultiplicative(ndemes, scaling,
final_multiplicative_trait(),
gvalue_to_fitness, noise);
return fwdpy11::multiplicative_trait_model(
ndemes, scaling, gvalue_to_fitness, noise);
}
return DiploidMultiplicative(ndemes, scaling,
final_multiplicative_fitness(),
gvalue_to_fitness, noise);
return fwdpy11::multiplicative_fitness_model(ndemes, scaling, noise);
}),
py::arg("scaling"), py::arg("gvalue_to_fitness"), py::arg("noise"),
py::arg("ndemes"));
Expand Down
78 changes: 6 additions & 72 deletions cpp_neutral_benchmark/cpp_neutral_benchmark.cc
Original file line number Diff line number Diff line change
Expand Up @@ -6,6 +6,7 @@
#include <fwdpy11/genetic_value_to_fitness/GeneticValueIsTrait.hpp>
#include <fwdpy11/genetic_values/dgvalue_pointer_vector.hpp>
#include <fwdpy11/genetic_values/fwdpp_wrappers/fwdpp_genetic_value.hpp>
#include <fwdpy11/genetic_values/DiploidMultiplicative.hpp>
#include <fwdpy11/regions/RecombinationRegions.hpp>
#include <fwdpy11/regions/MutationRegions.hpp>
#include <fwdpy11/evolvets/SampleRecorder.hpp>
Expand All @@ -21,73 +22,7 @@ namespace po = boost::program_options;
// the Python library C++ code. They are needed
// to define how Multiplicative works.

struct single_deme_multiplicative_het
{
inline void
operator()(double& d, const fwdpy11::Mutation& m) const
{
d *= (1. + m.s * m.h);
}
};

struct multi_deme_multiplicative_het
{
inline void
operator()(const std::size_t deme, double& d, const fwdpy11::Mutation& m) const
{
d *= (1. + m.esizes[deme] * m.heffects[deme]);
}
};

struct single_deme_multiplicative_hom
{
double scaling;
single_deme_multiplicative_hom(double s) : scaling(s)
{
}

inline void
operator()(double& d, const fwdpy11::Mutation& m) const
{
d *= (1. + scaling * m.s);
}
};

struct multi_deme_multiplicative_hom
{
double scaling;
multi_deme_multiplicative_hom(double s) : scaling(s)
{
}

inline void
operator()(const std::size_t deme, double& d, const fwdpy11::Mutation& m) const
{
d *= (1. + scaling * m.esizes[deme]);
}
};

struct final_multiplicative_trait
{
inline double
operator()(double d) const
{
return d - 1.0;
}
};

struct final_multiplicative_fitness
{
inline double
operator()(double d) const
{
return std::max(0.0, d);
}
};

using DiploidMultiplicative = fwdpy11::stateless_site_dependent_genetic_value_wrapper<
single_deme_multiplicative_het, single_deme_multiplicative_hom,
multi_deme_multiplicative_het, multi_deme_multiplicative_hom, 1>;
using fwdpy11::DiploidMultiplicative;

struct RecordNothing
// copied from internal code
Expand Down Expand Up @@ -175,8 +110,9 @@ simulate(const command_line_options& options)

fwdpy11_core::ForwardDemesGraph forward_demes_graph(o.str(), options.nsteps);

DiploidMultiplicative fitness(1, 2., final_multiplicative_fitness(), nullptr,
nullptr);
DiploidMultiplicative fitness(
1, 2., fwdpy11::final_multiplicative_fitness(),
[](const double) { return false; }, nullptr, nullptr);
fwdpy11::dgvalue_pointer_vector_ gvalue_pointers(fitness);

fwdpy11::no_ancient_samples no_ancient_samples{};
Expand All @@ -199,9 +135,7 @@ simulate(const command_line_options& options)

evolve_with_tree_sequences(
rng, pop, sample_recorder, options.simplification_interval, forward_demes_graph,
options.nsteps, 0., 0., mmodel, genetic_map, gvalue_pointers,
[](const fwdpy11::DiploidPopulation& /*pop*/, fwdpy11::SampleRecorder& /*sr*/) {
},
options.nsteps, 0., 0., mmodel, genetic_map, gvalue_pointers, no_ancient_samples,
stopping_criteron, record_nothing, tsoptions);
std::cout << pop.generation << ' ' << pop.N << ' ' << pop.tables->edges.size() << ' '
<< pop.tables->nodes.size() << '\n';
Expand Down
1 change: 1 addition & 0 deletions cpptests/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -8,6 +8,7 @@ SET(CPPTEST_SOURCES
forward_demes_graph_fixtures.cc
test_evolvets.cc
test_core_genetic_map_regions.cc
test_site_dependent_genetic_value.cc
)

add_executable(fwdpy11_cpp_tests ${CPPTEST_SOURCES})
Expand Down
Loading

0 comments on commit 0d36339

Please sign in to comment.