individual.hpp

Go to the documentation of this file.

#pragma once
#ifndef INDIVIDUAL_HPP_
#define INDIVIDUAL_HPP_

#include <cmath>
#include <iosfwd>
#include <string>
#include <vector>
#include <map>
#include <bitset>

namespace wtl {class sfmt19937_64;}

namespace edal {

namespace trait {
enum {
    toepad_size,
    limb_length,

    height_preference,
    diameter_preference,

    male_trait,
    female_trait,
    choosiness,

    neutral,

    size
};
} // namespace trait



struct IndividualParams {


    double BETA_PARAM = 3.0;

    size_t CARRYING_CAPACITY = 1000u;

    double AVG_NUM_OFFSPINRGS = 4.0;

    double HEIGHT_PREFERENCE = 2.0;

    double DIAMETER_PREFERENCE = 2.0;

    double TOEPAD_SELECTION = 2.0;

    double LIMB_SELECTION = 2.0;

    double PREF_COMPETITION = 2.0;

    double MORPH_COMPETITION = 2.0;

    double MATING_SIGMA = 0.05;


    double MU_LOCUS = 1e-4;


    unsigned long MUTATION_MASK = 0;

    double MIGRATION_RATE = 0.005;

    constexpr static double NORMAL_SIGMA = 0.3;

    constexpr static double C0 = 1.0;

    constexpr static double C1 = 0.5;
};

class Individual {
  public:
    using param_type = IndividualParams;
    using URBG = wtl::sfmt19937_64;

  private:
    constexpr static size_t NUM_LOCI_ = 8;

    constexpr static unsigned long FULL_BITS = (1u << NUM_LOCI_) - 1u;

    constexpr static unsigned long HALF_BITS = (1u << (NUM_LOCI_ / 2u)) - 1u;

    constexpr static double INV_NUM_LOCI_ = 0.5 / NUM_LOCI_;

    constexpr static size_t NUM_STEPS_ = 32;

    static std::map<std::vector<double>, double> KE_CACHE_;

  public:
    using Loci = std::bitset<NUM_LOCI_>;

    Individual(): Individual{
        {HALF_BITS, HALF_BITS, HALF_BITS, HALF_BITS, HALF_BITS, HALF_BITS, HALF_BITS, HALF_BITS},
        {HALF_BITS, HALF_BITS, HALF_BITS, HALF_BITS, HALF_BITS, HALF_BITS, HALF_BITS, HALF_BITS}} {}

    Individual(const std::vector<Loci>& egg,
               const std::vector<Loci>& sperm,
               const bool is_male=false):
        genotype_{egg, sperm},
        phenotype_(calc_phenotype()),
        sex_(is_male),
        ke_(effective_carrying_capacity_cache()) {}

    Individual(const std::vector<unsigned long>&);



    double effective_carrying_capacity_quad_unnormalized() const;


    double effective_carrying_capacity_exp_unnormalized() const;


    double effective_carrying_capacity_old_exp_unnormalized() const;

    double effective_carrying_capacity_cache() const;


    double preference_overlap(const Individual& other) const;


    double morphology_overlap(const Individual& other) const;


    double resource_overlap(const Individual& other) const {
        return std::exp(preference_overlap(other) + morphology_overlap(other));
    }


    double survival_probability(const double effective_num_competitors) const;


    double mating_probability(const Individual& male) const;


    double mating_probability_debarre(const Individual& male) const;


    double mating_probability_TPG2013(const Individual& male) const;


    std::vector<Loci> gametogenesis(URBG&) const;


    bool is_male() const {return sex_;}

    void change_sex() {sex_ = !sex_;}

    bool operator<(const Individual& other) const {
        return genotype_ < other.genotype_;
    }

    static std::string header();

    std::string str_detail() const;

    friend std::ostream& operator<< (std::ostream& ost, const Individual& ind);

    static void write_resource_abundance();
    static std::string possible_phenotypes();
    static std::string possible_geographic();

    static void param(const param_type& p) {PARAM_ = p;}
    static const param_type& param() {return PARAM_;}

  private:


    double habitat_preference_exp(const double height, const double diameter) const;


    double habitat_preference_quadratic(const double height, const double diameter) const;


    double calc_DI_analytical() const;


    double calc_DI_numerical() const;


    double calc_Dxi_analytical() const;


    double calc_Dxi_numerical() const;


    double fitness(const double height, const double diameter) const;



    std::vector<double> calc_phenotype() const {
        constexpr size_t n = trait::size;
        std::vector<double> output(n);
        for (size_t i=0; i<n; ++i) {
            output[i] = (genotype_.first[i].count() + genotype_.second[i].count()) * INV_NUM_LOCI_;
        }
        return output;
    }

    std::vector<double> intermediate_phenotypes() const;

    static const std::vector<std::string> INTERMEDIATE_KEYS_;

    friend double pdf_beta(const double height, const double diameter);
    friend double pdf_normal(const double height, const double diameter);
    friend double pdf_exp(const double height, const double diameter);
    template <class Func> friend double integrate_triangle(Func&& func);
    template <class Func> friend double integrate_square(Func&& func);

    // data member

    std::pair<std::vector<Loci>, std::vector<Loci> > genotype_;

    std::vector<double> phenotype_;

    bool sex_ = false;

    double ke_;

    static param_type PARAM_;
};

extern std::ostream& operator<<(std::ostream& ost, const Individual::Loci& bs);

} // namespace edal

namespace std {
inline bool operator< (const edal::Individual::Loci& lhs, const edal::Individual::Loci& rhs) {
    return lhs.to_ulong() < rhs.to_ulong();
}
}

#endif /* INDIVIDUAL_HPP_ */