Nemo  2.3.56
Simulate forward-in-time genetic evolution in a spatially explicit, individual-based stochastic simulator
LCE_Selection_base Class Reference

Base class performing (viability) selection on an arbitrary trait. More...

#include <LCEselection.h>

+ Inheritance diagram for LCE_Selection_base:
+ Collaboration diagram for LCE_Selection_base:

Public Member Functions

 LCE_Selection_base ()
 
virtual ~LCE_Selection_base ()
 
bool setSelectionMatrix ()
 
bool setLocalOptima ()
 
bool set_sel_model ()
 
bool set_fit_model ()
 
bool set_local_optima ()
 
bool set_param_rate_of_change ()
 
void set_std_rate_of_change ()
 
void addPhenotypicSD (unsigned int deme, double *stDev)
 
void changeLocalOptima ()
 
void resetCounters ()
 Resets the fitness counters. More...
 
void setMeans (unsigned int tot_ind)
 Computes the average fitness of each pedigree class. More...
 
double getMeanFitness (age_idx age)
 Computes the mean fitness of the whole population for a given age class. More...
 
double getMeanPatchFitness (age_idx age, unsigned int p)
 Computes the mean fitness in a given patch for a given age class. More...
 
double getMaxFitness (age_idx age)
 Computes the maximum fitness value of the whole population for a given age class. More...
 
double getMaxPatchFitness (age_idx age, unsigned int p)
 Computes the maximum fitness value in a given patch for a given age class. More...
 
double setMeanFitness (age_idx age)
 Sets the _mean_fitness variable to the value of the mean population fitness. More...
 
double getFitness (Individual *ind, unsigned int patch)
 Calls the fitness function according to the fitness model. More...
 
double getFitnessFixedEffect (Individual *offsprg, unsigned int patch, unsigned int trait)
 Returns the fitness of an individual in the fixed selection model. More...
 
double getFitnessDirect (Individual *offsprg, unsigned int patch, unsigned int trait)
 Returns the fitness of an individual following the direct selection model. More...
 
double getFitnessUnivariateQuadratic (Individual *ind, unsigned int patch, unsigned int trait)
 Quadratic fitness surface, approximates the Gaussian model for weak selection and/or small deviation from the optimum. More...
 
double getFitnessMultivariateGaussian (Individual *offsprg, unsigned int patch, unsigned int trait)
 Returns the fitness of an individual following the Gaussian selection model with one trait under selection. More...
 
double getFitnessUnivariateGaussian (Individual *offsprg, unsigned int patch, unsigned int trait)
 Returns the fitness of an individual following the Gaussian selection model with several traits under selection. More...
 
double getFitnessMultivariateGaussian_VE (Individual *offsprg, unsigned int patch, unsigned int trait)
 Returns the fitness of an individual following the Gaussian selection model with one trait under selection and environmental variance. More...
 
double getFitnessUnivariateGaussian_VE (Individual *offsprg, unsigned int patch, unsigned int trait)
 Returns the fitness of an individual following the Gaussian selection model with several traits under selection and environmental variance. More...
 
double getFitnessAbsolute (Individual *ind, unsigned int patch)
 Returns the raw fitness of the individual, without adjustment (absolute fitness). More...
 
double getFitnessRelative (Individual *ind, unsigned int patch)
 Returns the relative fitness of the individual, adjusted by a scaling factor. More...
 
void setScalingFactorLocal (age_idx age, unsigned int p)
 Sets the fitness scaling factor equal to the inverse of the mean local patch fitness. More...
 
void setScalingFactorGlobal (age_idx age, unsigned int p)
 Sets the fitness scaling factor equal to the inverse of the mean population fitness. More...
 
void setScalingFactorMaxLocal (age_idx age, unsigned int p)
 Sets the fitness scaling factor equal to the inverse of the maximum local patch fitness value. More...
 
void setScalingFactorMaxGlobal (age_idx age, unsigned int p)
 Sets the fitness scaling factor equal to the inverse of the maximum population fitness value. More...
 
void setScalingFactorAbsolute (age_idx age, unsigned int p)
 Resets the fitness scaling factor equal to one. More...
 
void doViabilitySelection (sex_t SEX, age_idx AGE, Patch *patch, unsigned int p)
 Selectively removes individuals in the population depending on their fitness. More...
 
Implementations
virtual bool setParameters ()
 
virtual void execute ()
 
virtual void loadStatServices (StatServices *loader)
 
virtual void loadFileServices (FileServices *loader)
 
virtual bool resetParameterFromSource (std::string param, SimComponent *cmpt)
 
virtual LifeCycleEventclone ()
 
virtual age_t removeAgeClass ()
 
virtual age_t addAgeClass ()
 
virtual age_t requiredAgeClass ()
 
- Public Member Functions inherited from LifeCycleEvent
 LifeCycleEvent (const char *name, const char *trait_link)
 Cstor. More...
 
virtual ~LifeCycleEvent ()
 
virtual void init (Metapop *popPtr)
 Sets the pointer to the current Metapop and the trait link if applicable. More...
 
virtual bool attach_trait (string trait)
 
virtual void set_paramset (std::string name, bool required, SimComponent *owner)
 
virtual void set_event_name (std::string &name)
 Set the name of the event (name of the ParamSet) and add the corresponding parameter to the set. More...
 
virtual void set_event_name (const char *name)
 
virtual string & get_event_name ()
 Accessor to the LCE's name. More...
 
virtual int get_rank ()
 Accessor to the LCE rank in the life cycle. More...
 
virtual void set_pop_ptr (Metapop *popPtr)
 Accessors for the population pointer. More...
 
virtual Metapopget_pop_ptr ()
 
- Public Member Functions inherited from SimComponent
 SimComponent ()
 
virtual ~SimComponent ()
 
virtual void loadUpdaters (UpdaterServices *loader)
 Loads the parameters and component updater onto the updater manager. More...
 
virtual void set_paramset (ParamSet *paramset)
 Sets the ParamSet member. More...
 
virtual void set_paramsetFromCopy (const ParamSet &PSet)
 Reset the set of parameters from a another set. More...
 
virtual ParamSetget_paramset ()
 ParamSet accessor. More...
 
virtual void add_parameter (Param *param)
 Interface to add a parameter to the set. More...
 
virtual void add_parameter (std::string Name, param_t Type, bool isRequired, bool isBounded, double low_bnd, double up_bnd)
 Interface to add a parameter to the set. More...
 
virtual void add_parameter (std::string Name, param_t Type, bool isRequired, bool isBounded, double low_bnd, double up_bnd, ParamUpdaterBase *updater)
 Interface to add a parameter and its updater to the set. More...
 
virtual Paramget_parameter (std::string name)
 Param getter. More...
 
virtual double get_parameter_value (std::string name)
 Param value getter. More...
 
virtual string get_name ()
 Returnd the name of the ParamSet, i.e. More...
 
virtual bool has_parameter (std::string name)
 Param getter. More...
 

Protected Attributes

double _letheq
 
double _base_fitness
 
double _mean_fitness
 
double _max_fitness
 
double _scaling_factor
 
bool _is_local
 
bool _is_absolute
 
double _fitness [5]
 Fitness counters, one for each pedigree class. More...
 
double _survival [5]
 
double _ind_cntr [5]
 
vector< string > _Traits
 The list of trait types under selection. More...
 
vector< unsigned int_TraitIndices
 The indices of the traits under selection. More...
 
vector< string > _SelectionModels
 The selection models associated with each trait under selection. More...
 
LCE_SelectionSH_stater
 
LCE_SelectionFH_writer
 
Quantitative traits parameters
int _selectTraitDimension
 Number of quantitative traits under selection. More...
 
vector< double > _selection_variance
 Patch-specific selection variance. More...
 
double _eVariance
 Evironmental variance. More...
 
Pointers to functions, set in LCE_Selection_base::setParameters
vector< double(LCE_Selection_base::*)(Individual *, unsigned int, unsigned int) > _getRawFitness
 A vector containing pointers to fitness function related to each trait under selection. More...
 
double(LCE_Selection_base::* _getFitness )(Individual *, unsigned int)
 Pointer to the function returning the individual fitness. More...
 
void(LCE_Selection_base::* _setScalingFactor )(age_idx, unsigned int)
 Pointer to the function used to set the fitness scaling factor when fitness is relative. More...
 
void(LCE_Selection_base::* _setNewLocalOptima )(void)
 Pointer to the function used to change the local phenotypic optima or its rate of change. More...
 
- Protected Attributes inherited from LifeCycleEvent
std::string _event_name
 The param name to be read in the init file. More...
 
Metapop_popPtr
 The ptr to the current Metapop. More...
 
std::string _LCELinkedTraitType
 The name of the linked trait. More...
 
int _LCELinkedTraitIndex
 The index in the individual's trait table of the linked trait. More...
 
- Protected Attributes inherited from SimComponent
ParamSet_paramSet
 The parameters container. More...
 

Private Attributes

Gaussian selection variables
vector< TMatrix * > _selection_matrix
 
vector< gsl_matrix * > _gsl_selection_matrix
 
gsl_vector * _diffs
 
gsl_vector * _res1
 
TMatrix_local_optima
 
TMatrix _rate_of_change_local_optima
 
bool _do_change_local_opt
 
bool _rate_of_change_is_std
 
unsigned int _set_std_rate_at_generation
 
int _std_rate_reference_patch
 
double * _phe
 Array to temporarily hold the phenotypic values of an individual, is never allocated. More...
 
Fixed selection (lethal-equivalents-based) parameters
double _Fpedigree [5]
 Array of pedigree values used in the fixed-selection model. More...
 
double _FitnessFixModel [5]
 Absolute fitness values of the five pedigree class for the fixed selection model (lethal equivalents model). More...
 

Friends

class LCE_SelectionSH
 The StatHandler associated class is a friend. More...
 
class LCE_SelectionFH
 
class LCE_Breed_Selection
 

Detailed Description

Base class performing (viability) selection on an arbitrary trait.

Implements several types of selection models, from fixed inbreeding depression-type model to Gaussian stabilizing selection on an adaptive landscape. This class is the base-class for the composite LCEs breed_selection and breed_selection_disperse.

Constructor & Destructor Documentation

◆ LCE_Selection_base()

LCE_Selection_base::LCE_Selection_base ( )
46 : LifeCycleEvent("viability_selection", ""),
50{
51 // mandatory parameter (no default):
52 add_parameter("selection_trait", STR, true, false, 0, 0); //no updaters here, to keep it safe...
53
56
57 add_parameter("selection_fitness_model",STR,false,false,0,0, updater);
58
60
61 add_parameter("selection_model",STR,true,false,0,0, updater); //mandatory, depends on trait
62 add_parameter("selection_matrix",MAT,false,false,0,0, updater);
63 add_parameter("selection_variance",DBL,false,false,0,0, updater);
64 add_parameter("selection_correlation",DBL,false,false,0,0, updater);
65 add_parameter("selection_trait_dimension",INT,false,false,0,0, updater);
66 add_parameter("selection_base_fitness",DBL,false,true,0,1, updater);
67 add_parameter("selection_lethal_equivalents",DBL,false,false,0,0, updater);
68 add_parameter("selection_pedigree_F",MAT,false,false,0,0, updater);
69 add_parameter("selection_randomize",BOOL,false,false,0,0, updater);
70 add_parameter("selection_environmental_variance", DBL, false, false, 0, 0, updater);
71
73 add_parameter("selection_local_optima",DBL,false,false,0,0, updater);
74
76 add_parameter("selection_rate_environmental_change", DBL, false, false, 0, 0, updater);
77 add_parameter("selection_std_rate_environmental_change", DBL, false, false, 0, 0, updater);
78 add_parameter("selection_std_rate_set_at_generation", INT, false, false, 0, 0, updater);
79 add_parameter("selection_std_rate_reference_patch", INT, false, false, 0, 0, updater);
80
81
82 add_parameter("selection_output", BOOL, false, false, 0, 0, 0);
83 add_parameter("selection_output_logtime", INT, false, false, 0, 0, 0);
84 add_parameter("selection_output_dir", STR, false, false, 0, 0, 0);
85}
int _selectTraitDimension
Number of quantitative traits under selection.
Definition: LCEselection.h:98
double _base_fitness
Definition: LCEselection.h:88
double(LCE_Selection_base::* _getFitness)(Individual *, unsigned int)
Pointer to the function returning the individual fitness.
Definition: LCEselection.h:123
double _max_fitness
Definition: LCEselection.h:89
vector< gsl_matrix * > _gsl_selection_matrix
Definition: LCEselection.h:56
LCE_SelectionFH * _writer
Definition: LCEselection.h:136
gsl_vector * _diffs
Definition: LCEselection.h:57
bool set_param_rate_of_change()
Definition: LCEselection.cc:605
bool _is_absolute
Definition: LCEselection.h:90
double _eVariance
Evironmental variance.
Definition: LCEselection.h:104
double _scaling_factor
Definition: LCEselection.h:89
bool _is_local
Definition: LCEselection.h:90
gsl_vector * _res1
Definition: LCEselection.h:57
bool set_sel_model()
Definition: LCEselection.cc:263
double _mean_fitness
Definition: LCEselection.h:89
bool set_fit_model()
Definition: LCEselection.cc:205
LCE_SelectionSH * _stater
Definition: LCEselection.h:135
double * _phe
Array to temporarily hold the phenotypic values of an individual, is never allocated.
Definition: LCEselection.h:64
vector< double(LCE_Selection_base::*)(Individual *, unsigned int, unsigned int) > _getRawFitness
A vector containing pointers to fitness function related to each trait under selection.
Definition: LCEselection.h:119
bool set_local_optima()
Definition: LCEselection.cc:575
TMatrix * _local_optima
Definition: LCEselection.h:58
vector< TMatrix * > _selection_matrix
Definition: LCEselection.h:55
Base class of the Life Cycle Events, declares the LCE interface.
Definition: lifecycleevent.h:73
Implementation of the ParamUpdaterBase interface.
Definition: param.h:363
virtual void add_parameter(Param *param)
Interface to add a parameter to the set.
Definition: simcomponent.h:112
@ DBL
Definition: types.h:78
@ MAT
Definition: types.h:78
@ BOOL
Definition: types.h:78
@ STR
Definition: types.h:78
@ INT
Definition: types.h:78

References SimComponent::add_parameter(), BOOL, DBL, INT, MAT, set_fit_model(), set_local_optima(), set_param_rate_of_change(), set_sel_model(), and STR.

Referenced by clone().

+ Here is the caller graph for this function:

◆ ~LCE_Selection_base()

LCE_Selection_base::~LCE_Selection_base ( )
virtual
90{
91 vector< TMatrix* >::iterator selIT = _selection_matrix.begin();
92 for(; selIT != _selection_matrix.end(); ++selIT)
93 if((*selIT)) delete (*selIT);
94 _selection_matrix.clear();
95
96 for(unsigned int i = 0; i < _gsl_selection_matrix.size(); ++i)
97 if(_gsl_selection_matrix[i]) gsl_matrix_free(_gsl_selection_matrix[i]);
99
100 if(_local_optima) delete _local_optima;
101 if(_diffs) gsl_vector_free(_diffs);
102 if(_res1) gsl_vector_free(_res1);
103 if(_stater) delete _stater;
104}

References _diffs, _gsl_selection_matrix, _local_optima, _res1, _selection_matrix, and _stater.

Member Function Documentation

◆ addAgeClass()

virtual age_t LCE_Selection_base::addAgeClass ( )
inlinevirtual

Implements LifeCycleEvent.

Reimplemented in LCE_Breed_Selection, and LCE_Breed_Selection_Disperse.

287{return NONE;}
#define NONE
No age flag.
Definition: types.h:48

References NONE.

◆ addPhenotypicSD()

void LCE_Selection_base::addPhenotypicSD ( unsigned int  deme,
double *  stDev 
)
1084{
1085 Patch *patch = _popPtr->getPatch(deme);
1086 Individual* ind;
1087 unsigned int table_size = patch->size(OFFSx);
1088
1089 double **phenot = new double* [_selectTraitDimension];
1090
1091 for (int i = 0; i < _selectTraitDimension; ++i) {
1092
1093 phenot[i] = new double [table_size];
1094
1095 }
1096
1097 unsigned int pos = 0;
1098
1099 for (unsigned int j = 0; j < patch->size(FEM, OFFSx) && pos < table_size; ++j) {
1100
1101 ind = patch->get(FEM, OFFSx, j);
1102 _phe = (double*)ind->getTraitValue(_LCELinkedTraitIndex);
1103
1104 for (int i = 0; i < _selectTraitDimension; ++i)
1105 phenot[i][pos] = _phe[i];
1106
1107 pos++;
1108
1109 }
1110
1111 for (unsigned int j = 0; j < patch->size(MAL, OFFSx) && pos < table_size; ++j) {
1112
1113 ind = patch->get(MAL, OFFSx, j);
1114 _phe = (double*)ind->getTraitValue(_LCELinkedTraitIndex);
1115
1116 for (int i = 0; i < _selectTraitDimension; ++i)
1117 phenot[i][pos] = _phe[i];
1118
1119 pos++;
1120
1121 }
1122
1123 assert(pos == table_size);
1124
1125 for (int i = 0; i < _selectTraitDimension; ++i) {
1126 stDev[i] += sqrt( my_variance_with_fixed_mean( phenot[i], table_size, my_mean(phenot[i], table_size) ) );
1127 delete [] phenot[i];
1128 }
1129
1130 delete [] phenot;
1131}
This class contains traits along with other individual information (sex, pedigree,...
Definition: individual.h:49
void * getTraitValue(IDX T)
Accessor to the value (phenotype) of a particular trait.
Definition: individual.h:271
int _LCELinkedTraitIndex
The index in the individual's trait table of the linked trait.
Definition: lifecycleevent.h:89
Metapop * _popPtr
The ptr to the current Metapop.
Definition: lifecycleevent.h:81
Patch * getPatch(unsigned int i)
Patch accessor, return the ith+1 patch in the metapop.
Definition: metapop.h:257
Second class in the metapopulation design structure, between the Metapop and Individual classes.
Definition: metapop.h:430
unsigned int size(age_t AGE)
Returns the size of the container of the appropriate age class(es) for both sexes.
Definition: metapop.h:496
Individual * get(sex_t SEX, age_idx AGE, unsigned int at)
Returns a pointer to the individual sitting at the index passed.
Definition: metapop.h:532
@ FEM
Definition: types.h:37
@ MAL
Definition: types.h:37
@ OFFSx
Definition: types.h:42
double my_mean(double *data, unsigned int size)
Definition: utils.cc:37
double my_variance_with_fixed_mean(double *data, unsigned int size, double mean)
Definition: utils.cc:62

References LifeCycleEvent::_LCELinkedTraitIndex, _phe, LifeCycleEvent::_popPtr, _selectTraitDimension, FEM, Patch::get(), Metapop::getPatch(), Individual::getTraitValue(), MAL, my_mean(), my_variance_with_fixed_mean(), OFFSx, and Patch::size().

Referenced by set_std_rate_of_change().

+ Here is the caller graph for this function:

◆ changeLocalOptima()

void LCE_Selection_base::changeLocalOptima ( )
1136{
1137 for (int i = 0; i < _selectTraitDimension; ++i)
1138 for (unsigned int j = 0; j < _local_optima->nrows(); ++j) {
1140 }
1141}
TMatrix _rate_of_change_local_optima
Definition: LCEselection.h:59
double get(unsigned int i, unsigned int j)
Accessor to element at row i and column j.
Definition: tmatrix.h:147
unsigned int nrows()
Definition: tmatrix.h:167
void plus(unsigned int i, unsigned int j, double value)
Adds a value to an element of the matrix.
Definition: tmatrix.h:210

References _local_optima, _rate_of_change_local_optima, _selectTraitDimension, TMatrix::get(), TMatrix::nrows(), and TMatrix::plus().

Referenced by set_param_rate_of_change(), and set_std_rate_of_change().

+ Here is the caller graph for this function:

◆ clone()

virtual LifeCycleEvent * LCE_Selection_base::clone ( )
inlinevirtual

Implements LifeCycleEvent.

Reimplemented in LCE_Breed_Selection, and LCE_Breed_Selection_Disperse.

285{return new LCE_Selection_base();}
LCE_Selection_base()
Definition: LCEselection.cc:46

References LCE_Selection_base().

◆ doViabilitySelection()

void LCE_Selection_base::doViabilitySelection ( sex_t  SEX,
age_idx  AGE,
Patch patch,
unsigned int  p 
)

Selectively removes individuals in the population depending on their fitness.

Calls the fitness function for each individual. Updates the fitness counters. The fitness scaling factor is set outside this function, in the execute() procedure.

Parameters
SEXthe sex class
AGEthe age-class index
patchthe focal patch
pthe focal patch index
940{
941 register Individual* ind;
942 register double fitness;
943 register unsigned int cat;
944
945 for(unsigned int i = 0; i < patch->size(SEX, AGE); i++) {
946
947 ind = patch->get(SEX, AGE, i);
948 cat = ind->getPedigreeClass();
949 fitness = getFitness( ind, p);
950
951 _fitness[cat] += fitness;
952 _ind_cntr[cat]++;
953
954 if(RAND::Uniform() > fitness ) {
955 //this individual dies
956 patch->remove(SEX, AGE, i);
957
958 _popPtr->recycle(ind);
959
960 i--;
961
962 } //else; this individual stays in the patch
963 else {
964 _survival[cat]++;
965 }
966
967 }
968}
void recycle(Individual *ind)
Put an individual in the recycling pool.
Definition: indfactory.h:62
unsigned int getPedigreeClass()
Returns the pedigree class of the individual, as set during offspring creation.
Definition: individual.h:179
double _fitness[5]
Fitness counters, one for each pedigree class.
Definition: LCEselection.h:92
double getFitness(Individual *ind, unsigned int patch)
Calls the fitness function according to the fitness model.
Definition: LCEselection.h:190
double _ind_cntr[5]
Definition: LCEselection.h:92
double _survival[5]
Definition: LCEselection.h:92
Individual * remove(sex_t SEX, age_idx AGE, unsigned int at)
Removes the individual sitting at the given index in the appropriate container.
Definition: metapop.h:573
static double Uniform()
Generates a random number from [0.0, 1.0[ uniformly distributed.
Definition: Uniform.h:125

References _fitness, _ind_cntr, LifeCycleEvent::_popPtr, _survival, Patch::get(), getFitness(), Individual::getPedigreeClass(), IndFactory::recycle(), Patch::remove(), Patch::size(), and RAND::Uniform().

Referenced by execute().

+ Here is the caller graph for this function:

◆ execute()

void LCE_Selection_base::execute ( )
virtual

Implements LifeCycleEvent.

Reimplemented in LCE_Breed_Selection, and LCE_Breed_Selection_Disperse.

898{
899 Patch * patch;
900 unsigned int popSize = _popPtr->size(OFFSPRG);
901
902#ifdef _DEBUG_
903 message("LCE_Selection_base::execute (tot offsprg nb: %i", popSize);
904#endif
905
907
909
910 if(_popPtr->getCurrentGeneration() == 1) {
911 set_local_optima(); //reset local optima to initial values
912 if(_rate_of_change_is_std) //reset rate of change relative to SD of that replicate
914 }
915
916 (this->*_setNewLocalOptima)();
917 }
918
919 for(unsigned int p = 0; p < _popPtr->getPatchNbr(); p++) {
920
921 (this->*_setScalingFactor)(OFFSx, p);
922
923 patch = _popPtr->getPatch(p);
924
925 doViabilitySelection(FEM, OFFSx, patch, p);
926
927 doViabilitySelection(MAL, OFFSx, patch, p);
928 }
929
930 setMeans(popSize);
931
932#ifdef _DEBUG_
933 message(", after selection: %i, mean fitness offspring=%f)\n",_popPtr->size(OFFSPRG), _mean_fitness);
934#endif
935}
void doViabilitySelection(sex_t SEX, age_idx AGE, Patch *patch, unsigned int p)
Selectively removes individuals in the population depending on their fitness.
Definition: LCEselection.cc:939
void(LCE_Selection_base::* _setScalingFactor)(age_idx, unsigned int)
Pointer to the function used to set the fitness scaling factor when fitness is relative.
Definition: LCEselection.h:128
bool _do_change_local_opt
Definition: LCEselection.h:60
void setMeans(unsigned int tot_ind)
Computes the average fitness of each pedigree class.
Definition: LCEselection.cc:983
void resetCounters()
Resets the fitness counters.
Definition: LCEselection.cc:972
void(LCE_Selection_base::* _setNewLocalOptima)(void)
Pointer to the function used to change the local phenotypic optima or its rate of change.
Definition: LCEselection.h:132
bool _rate_of_change_is_std
Definition: LCEselection.h:60
void set_std_rate_of_change()
Definition: LCEselection.cc:999
unsigned int size()
Get the total number of individuals present in the population, all sex and age classes together.
Definition: metapop.h:310
unsigned int getPatchNbr()
Definition: metapop.h:276
unsigned int getCurrentGeneration()
Definition: metapop.h:294
void message(const char *message,...)
Definition: output.cc:40
#define OFFSPRG
Offspring age class flag.
Definition: types.h:50

References _do_change_local_opt, _mean_fitness, LifeCycleEvent::_popPtr, _rate_of_change_is_std, _setNewLocalOptima, _setScalingFactor, doViabilitySelection(), FEM, Metapop::getCurrentGeneration(), Metapop::getPatch(), Metapop::getPatchNbr(), MAL, message(), OFFSPRG, OFFSx, resetCounters(), set_local_optima(), set_std_rate_of_change(), setMeans(), and Metapop::size().

◆ getFitness()

double LCE_Selection_base::getFitness ( Individual ind,
unsigned int  patch 
)
inline

Calls the fitness function according to the fitness model.

The fitness model can be "absolute", "relative_local" or "relative_global".

Parameters
indthe focal indvidual, we want to know its fitness
patchthe index of the patch of the focal individual
191 {
192 return (this->*_getFitness)(ind, patch);
193 }

References _getFitness.

Referenced by LCE_Breed_Selection::do_breed_selection_FecFitness(), doViabilitySelection(), getMaxPatchFitness(), getMeanFitness(), getMeanPatchFitness(), LCE_Breed_Selection::makeOffspringWithSelection(), LCE_SelectionSH::setDataTable(), and LCE_Breed_Selection::setReproScaledFitness_sum().

+ Here is the caller graph for this function:

◆ getFitnessAbsolute()

double LCE_Selection_base::getFitnessAbsolute ( Individual ind,
unsigned int  patch 
)

Returns the raw fitness of the individual, without adjustment (absolute fitness).

Calls the fitness function according to the right selection model.

Parameters
indthe focal indvidual, we want to know its fitness
patchthe index of the patch of the focal individual
780{
781 double fitness = 1;
782 //at this point, _getRawFitness.size() == _TraitIndices.size()
783 //and we assume that the functions are "aligned" with the traits
784 for(unsigned int i = 0; i < _getRawFitness.size(); i++)
785 fitness *= (this->*_getRawFitness[i])(ind, patch, _TraitIndices[i]);
786 return fitness;
787}
vector< unsigned int > _TraitIndices
The indices of the traits under selection.
Definition: LCEselection.h:111

References _getRawFitness, and _TraitIndices.

Referenced by getFitnessRelative(), and set_fit_model().

+ Here is the caller graph for this function:

◆ getFitnessDirect()

double LCE_Selection_base::getFitnessDirect ( Individual offsprg,
unsigned int  patch,
unsigned int  trait 
)
inline

Returns the fitness of an individual following the direct selection model.

203 {
204 return *(double*)offsprg->getTraitValue(trait);
205 }

References Individual::getTraitValue().

Referenced by set_sel_model().

+ Here is the caller graph for this function:

◆ getFitnessFixedEffect()

double LCE_Selection_base::getFitnessFixedEffect ( Individual offsprg,
unsigned int  patch,
unsigned int  trait 
)
inline

Returns the fitness of an individual in the fixed selection model.

197 {
198 return _FitnessFixModel[ offsprg->getPedigreeClass() ];
199 }
double _FitnessFixModel[5]
Absolute fitness values of the five pedigree class for the fixed selection model (lethal equivalents ...
Definition: LCEselection.h:76

References _FitnessFixModel, and Individual::getPedigreeClass().

Referenced by set_sel_model().

+ Here is the caller graph for this function:

◆ getFitnessMultivariateGaussian()

double LCE_Selection_base::getFitnessMultivariateGaussian ( Individual offsprg,
unsigned int  patch,
unsigned int  trait 
)

Returns the fitness of an individual following the Gaussian selection model with one trait under selection.

721{
722 register double res2;
723
724// if(!ind) fatal("passing NULL ind ptr to LCE_Selection_base::getFitnessMultivariateGaussian!!!\n");
725
726 _phe = (double*)ind->getTraitValue(trait);
727
728 for( int i = 0; i < _selectTraitDimension; i++)
729 gsl_vector_set(_diffs, i, _phe[i] - _local_optima->get(patch, i));
730
731 //(diff)T * W * diff:
732 //right partial product:
733 gsl_blas_dsymv(CblasUpper, 1.0, _gsl_selection_matrix[patch], _diffs, 0.0, _res1);
734 //left product:
735 gsl_blas_ddot(_diffs, _res1, &res2);
736
737 return exp( -0.5 * res2 );
738}

References _diffs, _gsl_selection_matrix, _local_optima, _phe, _res1, _selectTraitDimension, TMatrix::get(), and Individual::getTraitValue().

Referenced by set_sel_model().

+ Here is the caller graph for this function:

◆ getFitnessMultivariateGaussian_VE()

double LCE_Selection_base::getFitnessMultivariateGaussian_VE ( Individual offsprg,
unsigned int  patch,
unsigned int  trait 
)

Returns the fitness of an individual following the Gaussian selection model with one trait under selection and environmental variance.

759{
760 register double res2;
761
762 _phe = (double*)ind->getTraitValue(trait);
763
764 //add the environmental variance here:
765 for( int i = 0; i < _selectTraitDimension; i++)
766 gsl_vector_set(_diffs, i, _phe[i] + RAND::Gaussian(_eVariance) - _local_optima->get(patch, i));
767
768 //(diff)T * W * diff:
769 //right partial product:
770 gsl_blas_dsymv(CblasUpper, 1.0, _gsl_selection_matrix[patch], _diffs, 0.0, _res1);
771 //left product:
772 gsl_blas_ddot(_diffs, _res1, &res2);
773
774 return exp( -0.5 * res2 );
775}
static double Gaussian(double sigma)
Definition: Uniform.h:262

References _diffs, _eVariance, _gsl_selection_matrix, _local_optima, _phe, _res1, _selectTraitDimension, RAND::Gaussian(), TMatrix::get(), and Individual::getTraitValue().

Referenced by set_sel_model().

+ Here is the caller graph for this function:

◆ getFitnessRelative()

double LCE_Selection_base::getFitnessRelative ( Individual ind,
unsigned int  patch 
)
inline

Returns the relative fitness of the individual, adjusted by a scaling factor.

The scaling factor is set according to the relative fitness model (local or global) outside of this function. If the scaling factor is one, this is the absolute fitness of the individual. Calls the fitness function according to the right selection model.

Parameters
indthe focal indvidual, we want to know its fitness
patchthe index of the patch of the focal individual
236 {
237 return getFitnessAbsolute(ind, patch) * _scaling_factor;
238 }
double getFitnessAbsolute(Individual *ind, unsigned int patch)
Returns the raw fitness of the individual, without adjustment (absolute fitness).
Definition: LCEselection.cc:779

References _scaling_factor, and getFitnessAbsolute().

Referenced by set_fit_model().

+ Here is the caller graph for this function:

◆ getFitnessUnivariateGaussian()

double LCE_Selection_base::getFitnessUnivariateGaussian ( Individual offsprg,
unsigned int  patch,
unsigned int  trait 
)

Returns the fitness of an individual following the Gaussian selection model with several traits under selection.

706{
707 register double res2, diff;
708
709 _phe = (double*)ind->getTraitValue(trait);
710
711 diff = _phe[0] - _local_optima->get(patch, 0);
712
713 res2 = diff*diff / _selection_variance[patch];
714
715 return exp( -0.5 * res2 );
716}
vector< double > _selection_variance
Patch-specific selection variance.
Definition: LCEselection.h:101

References _local_optima, _phe, _selection_variance, TMatrix::get(), and Individual::getTraitValue().

Referenced by set_sel_model().

+ Here is the caller graph for this function:

◆ getFitnessUnivariateGaussian_VE()

double LCE_Selection_base::getFitnessUnivariateGaussian_VE ( Individual offsprg,
unsigned int  patch,
unsigned int  trait 
)

Returns the fitness of an individual following the Gaussian selection model with several traits under selection and environmental variance.

743{
744 register double res2, diff;
745
746 _phe = (double*)ind->getTraitValue(trait);
747
748 //add the environmental variance here:
749 diff = _phe[0] + RAND::Gaussian(_eVariance) - _local_optima->get(patch, 0);
750
751 res2 = diff*diff / _selection_variance[patch];
752
753 return exp( -0.5 * res2 );
754}

References _eVariance, _local_optima, _phe, _selection_variance, RAND::Gaussian(), TMatrix::get(), and Individual::getTraitValue().

Referenced by set_sel_model().

+ Here is the caller graph for this function:

◆ getFitnessUnivariateQuadratic()

double LCE_Selection_base::getFitnessUnivariateQuadratic ( Individual ind,
unsigned int  patch,
unsigned int  trait 
)

Quadratic fitness surface, approximates the Gaussian model for weak selection and/or small deviation from the optimum.

691{
692 register double res2, diff;
693
694 _phe = (double*)ind->getTraitValue(trait);
695
696 diff = _phe[0] - _local_optima->get(patch, 0);
697
698 res2 = diff*diff / _selection_variance[patch];
699
700 return 1 - res2;
701}

References _local_optima, _phe, _selection_variance, TMatrix::get(), and Individual::getTraitValue().

Referenced by set_sel_model().

+ Here is the caller graph for this function:

◆ getMaxFitness()

double LCE_Selection_base::getMaxFitness ( age_idx  age)

Computes the maximum fitness value of the whole population for a given age class.

Parameters
agethe age-class index
826{
827 double max = 0, local_max;
828
829 for(unsigned int i = 0, npatch = _popPtr->getPatchNbr(); i < npatch; i++) {
830
831 local_max = getMaxPatchFitness(age, i);
832
833 max = (local_max > max ? local_max : max);
834 }
835
836 return max;
837}
double getMaxPatchFitness(age_idx age, unsigned int p)
Computes the maximum fitness value in a given patch for a given age class.
Definition: LCEselection.cc:841

References LifeCycleEvent::_popPtr, getMaxPatchFitness(), and Metapop::getPatchNbr().

Referenced by setScalingFactorMaxGlobal().

+ Here is the caller graph for this function:

◆ getMaxPatchFitness()

double LCE_Selection_base::getMaxPatchFitness ( age_idx  age,
unsigned int  p 
)

Computes the maximum fitness value in a given patch for a given age class.

Parameters
agethe age-class index
pthe patch index
842{
843 double max = 0, fit;
844 Patch *patch = _popPtr->getPatch(p);
845
846 for(unsigned int j = 0, size = patch->size(FEM, age); j < size; j++) {
847 fit = getFitness( patch->get(FEM, age, j), p);
848 max = (fit > max ? fit : max);
849 }
850 for(unsigned int j = 0, size = patch->size(MAL, age); j < size; j++){
851 fit = getFitness( patch->get(MAL, age, j), p);
852 max = (fit > max ? fit : max);
853 }
854
855 return max;
856}

References LifeCycleEvent::_popPtr, FEM, Patch::get(), getFitness(), Metapop::getPatch(), MAL, and Patch::size().

Referenced by getMaxFitness(), and setScalingFactorMaxLocal().

+ Here is the caller graph for this function:

◆ getMeanFitness()

double LCE_Selection_base::getMeanFitness ( age_idx  age)

Computes the mean fitness of the whole population for a given age class.

Parameters
agethe age-class index
792{
793 double mean = 0;
794 Patch *patch;
795// age_idx age = (AGE == ADULTS ? ADLTx : OFFSx);
796
797 for(unsigned int i = 0, npatch = _popPtr->getPatchNbr(); i < npatch; i++) {
798 patch = _popPtr->getPatch(i);
799 for(unsigned int j = 0, size = patch->size(FEM, age); j < size; j++)
800 mean += getFitness( patch->get(FEM, age, j), i);
801 for(unsigned int j = 0, size = patch->size(MAL, age); j < size; j++)
802 mean += getFitness( patch->get(MAL, age, j), i);
803 }
804 return mean/_popPtr->size(age);
805}

References LifeCycleEvent::_popPtr, FEM, Patch::get(), getFitness(), Metapop::getPatch(), Metapop::getPatchNbr(), MAL, Metapop::size(), and Patch::size().

Referenced by LCE_Breed_Selection_Disperse::execute(), setMeanFitness(), and setScalingFactorGlobal().

+ Here is the caller graph for this function:

◆ getMeanPatchFitness()

double LCE_Selection_base::getMeanPatchFitness ( age_idx  age,
unsigned int  p 
)

Computes the mean fitness in a given patch for a given age class.

Parameters
agethe age-class index
pthe patch index
810{
811 double mean = 0;
812 Patch *patch = _popPtr->getPatch(p);
813
814 for(unsigned int j = 0, size = patch->size(FEM, age); j < size; j++)
815 mean += getFitness( patch->get(FEM, age, j), p);
816
817 for(unsigned int j = 0, size = patch->size(MAL, age); j < size; j++)
818 mean += getFitness( patch->get(MAL, age, j), p);
819
820 return mean/patch->size(age);
821}

References LifeCycleEvent::_popPtr, FEM, Patch::get(), getFitness(), Metapop::getPatch(), MAL, and Patch::size().

Referenced by LCE_Breed_Selection_Disperse::breed_selection_disperse(), and setScalingFactorLocal().

+ Here is the caller graph for this function:

◆ loadFileServices()

void LCE_Selection_base::loadFileServices ( FileServices loader)
virtual

Implements SimComponent.

Reimplemented in LCE_Breed_Selection, and LCE_Breed_Selection_Disperse.

117{
118 if(_paramSet->isSet("selection_output")) {
119
120 if(_writer == NULL) _writer = new LCE_SelectionFH(this);
121
122 Param* param = get_parameter("selection_output_logtime");
123
124 if(!param->isSet()) fatal("parameter \"selection_output_logtime\" is missing\n");
125
126 if(param->isMatrix()) {
127
128 TMatrix temp;
129
130 param->getMatrix(&temp);
131
132 _writer->set_multi(true, true, 1, &temp, get_parameter("selection_output_dir")->getArg());
133 // rpl_per, gen_per, rpl_occ, gen_occ, rank, path, self-ref
134 } else _writer->set(true, param->isSet(), 1, (param->isSet() ? (int)param->getValue() : 0), 0, get_parameter("selection_output_dir")->getArg(), this);
135
136
137 loader->attach(_writer);
138
139 } else {
140
141 if(_writer) delete _writer;
142
143 return;
144 }
145
146
147}
virtual void set(bool rpl_per, bool gen_per, int rpl_occ, int gen_occ, int rank, string path, LCE *event)
Definition: filehandler.h:265
virtual void set_multi(bool rpl_per, bool gen_per, int rpl_occ, TMatrix *Occ, string path)
Definition: filehandler.h:197
virtual void attach(Handler *FH)
Attaches the FileHandler to the current list (_writers) of the FileServices.
Definition: fileservices.cc:60
friend class LCE_SelectionFH
Definition: LCEselection.h:83
bool isSet()
Accessor to the status flag.
Definition: param.h:288
This structure stores one parameter, its definition and its string argument.
Definition: param.h:54
double getValue()
Returns the argument value according to its type.
Definition: param.cc:347
bool isMatrix()
Checks if the argument is of matrix type.
Definition: param.h:172
void getMatrix(TMatrix *mat)
Sets the matrix from the argument string if the parameter is set and of matrix type.
Definition: param.cc:357
bool isSet()
Definition: param.h:140
virtual Param * get_parameter(std::string name)
Param getter.
Definition: simcomponent.h:139
ParamSet * _paramSet
The parameters container.
Definition: simcomponent.h:48
A class to handle matrix in params, coerces matrix into a vector of same total size.
Definition: tmatrix.h:49
void fatal(const char *str,...)
Definition: output.cc:96

References SimComponent::_paramSet, _writer, FileServices::attach(), fatal(), SimComponent::get_parameter(), Param::getMatrix(), Param::getValue(), Param::isMatrix(), Param::isSet(), ParamSet::isSet(), LCE_SelectionFH, EventFileHandler< LCE >::set(), and FileHandler::set_multi().

◆ loadStatServices()

void LCE_Selection_base::loadStatServices ( StatServices loader)
virtual

Implements SimComponent.

Reimplemented in LCE_Breed_Selection, and LCE_Breed_Selection_Disperse.

109{
110 if(_stater == NULL) _stater = new LCE_SelectionSH(this);
111 loader->attach(_stater);
112}
friend class LCE_SelectionSH
The StatHandler associated class is a friend.
Definition: LCEselection.h:82
virtual void attach(Handler *H)
attach the StatHandler to the current list (_statHandlers) of the StatServices
Definition: statservices.cc:177

References _stater, StatServices::attach(), and LCE_SelectionSH.

Referenced by LCE_Breed_Selection::loadStatServices(), and LCE_Breed_Selection_Disperse::loadStatServices().

+ Here is the caller graph for this function:

◆ removeAgeClass()

virtual age_t LCE_Selection_base::removeAgeClass ( )
inlinevirtual

Implements LifeCycleEvent.

Reimplemented in LCE_Breed_Selection, and LCE_Breed_Selection_Disperse.

286{return NONE;}

References NONE.

◆ requiredAgeClass()

virtual age_t LCE_Selection_base::requiredAgeClass ( )
inlinevirtual

Implements LifeCycleEvent.

Reimplemented in LCE_Breed_Selection, and LCE_Breed_Selection_Disperse.

288{return OFFSPRG;}

References OFFSPRG.

◆ resetCounters()

void LCE_Selection_base::resetCounters ( )

Resets the fitness counters.

973{
974 for(unsigned int i = 0; i < 5; i++) {
975 _fitness[i] = 0;
976 _survival[i] = 0;
977 _ind_cntr[i] = 0;
978 }
979}

References _fitness, _ind_cntr, and _survival.

Referenced by LCE_Breed_Selection::execute(), LCE_Breed_Selection_Disperse::execute(), execute(), and setParameters().

+ Here is the caller graph for this function:

◆ resetParameterFromSource()

virtual bool LCE_Selection_base::resetParameterFromSource ( std::string  param,
SimComponent cmpt 
)
inlinevirtual

Implements SimComponent.

Reimplemented in LCE_Breed_Selection, and LCE_Breed_Selection_Disperse.

284{return false;}

◆ set_fit_model()

bool LCE_Selection_base::set_fit_model ( )
206{
207
208 if(_paramSet->isSet("selection_fitness_model")) {
209
210 string fit_model = _paramSet->getArg("selection_fitness_model");
211
212 if( fit_model == "absolute") {
213
214 _is_local = false;
215 _is_absolute = true;
218
219 } else if( fit_model == "relative_local" ) {
220
221 _is_local = true;
222 _is_absolute = false;
225
226 } else if(fit_model == "relative_global") {
227
228 _is_local = false;
229 _is_absolute = false;
232
233 } else if(fit_model == "relative_max_local") {
234
235 _is_local = true;
236 _is_absolute = false;
239
240 } else if(fit_model == "relative_max_global") {
241
242 _is_local = false;
243 _is_absolute = false;
246
247 } else {
248 return error("Unknown fitness model \"%s\"", fit_model.c_str());
249 }
250 } //default case:
251 else {
252 _is_local = false;
253 _is_absolute = true;
256 }
257
258 return true;
259}
double getFitnessRelative(Individual *ind, unsigned int patch)
Returns the relative fitness of the individual, adjusted by a scaling factor.
Definition: LCEselection.h:235
void setScalingFactorAbsolute(age_idx age, unsigned int p)
Resets the fitness scaling factor equal to one.
Definition: LCEselection.h:267
void setScalingFactorMaxGlobal(age_idx age, unsigned int p)
Sets the fitness scaling factor equal to the inverse of the maximum population fitness value.
Definition: LCEselection.cc:887
void setScalingFactorMaxLocal(age_idx age, unsigned int p)
Sets the fitness scaling factor equal to the inverse of the maximum local patch fitness value.
Definition: LCEselection.cc:879
void setScalingFactorGlobal(age_idx age, unsigned int p)
Sets the fitness scaling factor equal to the inverse of the mean population fitness.
Definition: LCEselection.cc:869
void setScalingFactorLocal(age_idx age, unsigned int p)
Sets the fitness scaling factor equal to the inverse of the mean local patch fitness.
Definition: LCEselection.cc:861
string getArg(string name)
Accessor to the parameters argument string.
Definition: param.h:300
int error(const char *str,...)
Definition: output.cc:77

References _getFitness, _is_absolute, _is_local, SimComponent::_paramSet, _setScalingFactor, error(), ParamSet::getArg(), getFitnessAbsolute(), getFitnessRelative(), ParamSet::isSet(), setScalingFactorAbsolute(), setScalingFactorGlobal(), setScalingFactorLocal(), setScalingFactorMaxGlobal(), and setScalingFactorMaxLocal().

Referenced by LCE_Selection_base(), and setParameters().

+ Here is the caller graph for this function:

◆ set_local_optima()

bool LCE_Selection_base::set_local_optima ( )
576{//set the traits' local optima, _selectTraitDimension must be set before that (= nbr of traits to select on)
577 string model = _paramSet->getArg("selection_model");
578
579 if( model == "fix" || model == "direct") return true;
580
581 if( (model == "gaussian" || model == "quadratic")
582 && !get_parameter("selection_local_optima")->isSet()) {
583
584 return error("parameter \"selection_local_optima\" must be set to have Gaussian or quadratic selection.\n");
585
586 } else {
587
588 TMatrix tmp_mat;
589
590 if(_local_optima == 0) _local_optima = new TMatrix();
591
593
594 _paramSet->getMatrix("selection_local_optima", &tmp_mat);
595
596 return setSpatialMatrix("selection_local_optima", "\"selection_trait_dimension\"", &tmp_mat, _local_optima,
597 _selectTraitDimension, _popPtr->getPatchNbr(), _paramSet->isSet("selection_randomize"));
598 }
599
600 return true;
601}
void getMatrix(string name, TMatrix *mat)
Accessor to the parameters matrix.
Definition: param.h:304
void reset(unsigned int rows, unsigned int cols)
Re-allocate the existing matrix with assigned rows and cols dimensions.
Definition: tmatrix.h:116
bool setSpatialMatrix(string param, string numColCondition, TMatrix *inMat, TMatrix *outMat, unsigned int nVal, unsigned int patchNbr, bool doRandomize)
Definition: utils.cc:114

References _local_optima, SimComponent::_paramSet, LifeCycleEvent::_popPtr, _selectTraitDimension, error(), SimComponent::get_parameter(), ParamSet::getArg(), ParamSet::getMatrix(), Metapop::getPatchNbr(), ParamSet::isSet(), TMatrix::reset(), and setSpatialMatrix().

Referenced by execute(), LCE_Selection_base(), and setParameters().

+ Here is the caller graph for this function:

◆ set_param_rate_of_change()

bool LCE_Selection_base::set_param_rate_of_change ( )
606{
608 _do_change_local_opt = false;
611
612
613 if (!get_parameter("selection_rate_environmental_change")->isSet()
614 && !get_parameter("selection_std_rate_environmental_change")->isSet() )
615 return true;
616
617 if (get_parameter("selection_rate_environmental_change")->isSet()
618 && get_parameter("selection_std_rate_environmental_change")->isSet() ) {
619 return error("both \"selection_rate_environmental_change\" and \"selection_std_rate_environmental_change\" are set, need only one.\n");
620 }
621
622 TMatrix tmpMat;
623
624 if (get_parameter("selection_rate_environmental_change")->isSet()) {
625
626 if(!get_parameter("selection_rate_environmental_change")->isMatrix()) {
627
628 double val = get_parameter_value("selection_rate_environmental_change");
629
630 tmpMat.reset(1, _selectTraitDimension);
631 tmpMat.assign(val);
632
633 } else {
634 get_parameter("selection_rate_environmental_change")->getMatrix(&tmpMat);
635 }
636
638
639 } else if (get_parameter("selection_std_rate_environmental_change")->isSet()){
640
641 if(!get_parameter("selection_std_rate_environmental_change")->isMatrix()) {
642
643 double val = get_parameter_value("selection_std_rate_environmental_change");
644
645 tmpMat.reset(1, _selectTraitDimension);
646 tmpMat.assign(val);
647
648 } else {
649 get_parameter("selection_std_rate_environmental_change")->getMatrix(&tmpMat);
650 }
651
653
654 if(get_parameter("selection_std_rate_set_at_generation")->isSet())
655 _set_std_rate_at_generation = (unsigned int)get_parameter_value("selection_std_rate_set_at_generation");
656 else
658
659 //check if phenotypic SD is to be computed in a single reference patch
660 //is -1 if parameter not set, which corresponds to the whole population then
661 _std_rate_reference_patch = get_parameter_value("selection_std_rate_reference_patch");
662
664 }
665
666 if(tmpMat.nrows() > _popPtr->getPatchNbr())
667 return error("The matrix of rate of change in local optima has more rows than the number of patches\n");
668
669 if((int)tmpMat.ncols() > _selectTraitDimension)
670 return error("The matrix of rate of change in local optima has more columns than the number of quantitative traits\n");
671
673
675
676 unsigned int nCol = tmpMat.ncols(), nRow = tmpMat.nrows();
677
678 //copy values, with pattern propagation
679 for (int p = 0; p < _popPtr->getPatchNbr(); ++p) {
680 for (int i = 0; i < _selectTraitDimension; ++i) {
681 _rate_of_change_local_optima.set(p, i, tmpMat.get(p % nRow, i % nCol));
682 }
683 }
684
685 return true;
686}
unsigned int _set_std_rate_at_generation
Definition: LCEselection.h:61
int _std_rate_reference_patch
Definition: LCEselection.h:62
void changeLocalOptima()
Definition: LCEselection.cc:1135
virtual double get_parameter_value(std::string name)
Param value getter.
Definition: simcomponent.h:143
unsigned int ncols()
Definition: tmatrix.h:170
void set(unsigned int i, unsigned int j, double val)
Sets element at row i and column j to value val.
Definition: tmatrix.h:102
void assign(double val)
Assigns a value to all element of the matrix.
Definition: tmatrix.h:110

References _do_change_local_opt, LifeCycleEvent::_popPtr, _rate_of_change_is_std, _rate_of_change_local_optima, _selectTraitDimension, _set_std_rate_at_generation, _setNewLocalOptima, _std_rate_reference_patch, TMatrix::assign(), changeLocalOptima(), error(), TMatrix::get(), SimComponent::get_parameter(), SimComponent::get_parameter_value(), Param::getMatrix(), Metapop::getPatchNbr(), TMatrix::ncols(), TMatrix::nrows(), TMatrix::reset(), TMatrix::set(), and set_std_rate_of_change().

Referenced by LCE_Selection_base(), and setParameters().

+ Here is the caller graph for this function:

◆ set_sel_model()

bool LCE_Selection_base::set_sel_model ( )
264{
265
266 _selectTraitDimension = (int)get_parameter_value("selection_trait_dimension");
267
268 _base_fitness = get_parameter_value("selection_base_fitness");
269
270 // -------------------------------------------------------------------------------------
271 // ENVIRONMENTAL VARIANCE
272
273 if(get_parameter("selection_environmental_variance")->isSet())
274 //the variable actually holds the standard dev...
275 _eVariance = sqrt(get_parameter_value("selection_environmental_variance"));
276 else
277 _eVariance = 0;
278
279 // -------------------------------------------------------------------------------------
280 // SELECTION MODELS
281
282 // empty containers:
283 _SelectionModels.clear();
284 _getRawFitness.clear();
285
286 // default case: direct fitness model
287 if(!_paramSet->isSet("selection_model")) {
288
290
291 _SelectionModels.push_back("direct");
292
293 } else {
294
295 _SelectionModels = get_parameter("selection_model")->getMultiArgs();
296 } //end_if isSet()
297
298 // check consistency with selection_trait info
299 if (_SelectionModels.size() != _TraitIndices.size())
300 return error("\"selection_trait\" and \"selection_model\" must have the same number of arguments.");
301
302 // -------------------------------------------------------------------------------------
303 // SET FUNCTION POINTERS TO FITNESS FUNCTIONS
304
305 string sel_model;
306
307 for(unsigned int t = 0; t < _SelectionModels.size(); t++) {
308
309 sel_model = _SelectionModels[t];
310
311
312 // FIX -----------------------------------------------------------------------
313 if(sel_model == "fix") {
314
315 //check if the trait type in the trait table is correct
316 if(_Traits[t] != "fix")
317 return error("the selection model \"fix\" does not match with the trait type \"%s\"\n", _Traits[t].c_str());
318
319 if(!get_parameter("selection_lethal_equivalents")->isSet()) {
320
321 return error("\"selection_lethal_equivalents\" parameter is missing with \"fix\" selection!\n");
322
323 } else
324 _letheq = get_parameter_value("selection_lethal_equivalents");
325
326 if(!get_parameter("selection_base_fitness")->isSet()) {
327
328 warning("\"selection_base_fitness\" parameter is missing under fix selection model, setting it to 1.\n");
329
330 _base_fitness = 1.0;
331 }
332
333 if(!get_parameter("selection_pedigree_F")->isSet()) {
334
335 return error("\"selection_pedigree_F\" parameter is missing with \"fix\" selection!\n");
336
337 } else {
338
339 TMatrix tmp_mat;
340
341 get_parameter("selection_pedigree_F")->getMatrix(&tmp_mat);
342
343 if(tmp_mat.getNbCols() != 5) {
344 return error("\"selection_pedigree_F\" must be an array of size 5.\n");
345 }
346
347 for(unsigned int i = 0; i < 5; i++) _Fpedigree[i] = tmp_mat.get(0,i);
348 }
349
350 for(unsigned int i = 0; i < 5; i++)
352
353 //everything has been set correctly, add function pointer to container
355
356 // DIRECT -----------------------------------------------------------------------
357 } else if(sel_model == "direct") {
358
359 //check if the trait type in the trait table is correct
360 if(_Traits[t] != "delet" && _Traits[t] != "dmi")
361 return error("the selection model \"direct\" does not match with the trait type \"%s\"\n", _Traits[t].c_str());
362
364
365 // QUADRATIC --------------------------------------------------------------------
366 } else if(sel_model == "quadratic") {
367
368 //check if the trait type in the trait table is correct
369 if(_Traits[t] != "quant")
370 return error("the selection model \"quadratic\" does not match with the trait type \"%s\"\n", _Traits[t].c_str());
371
372 if( get_parameter("selection_trait_dimension")->isSet())
373 _selectTraitDimension = get_parameter_value("selection_trait_dimension");
374 else
376
377 if(_selectTraitDimension == 1){
378
379 //check if the number of quanti traits modeled match with the dimensionality:
380 if((int)_popPtr->getTraitPrototype("quant")->get_parameter_value("quanti_traits") != _selectTraitDimension)
381 return error("The number of quantitative traits in \"quanti_traits\" does not match with\
382 the selection trait dimensionality in \"selection_trait_dimension\"\n");
383
384 if(!setSelectionMatrix()) return false; //this to set the selection variance params
385
386 //now add the pointer to fitness function:
388
389 } else return error("\"quadratic\" fitness model implemented for a single trait only.\n");
390
391
392 // GAUSSIAN --------------------------------------------------------------------
393 } else if(sel_model == "gaussian") {
394
395 //check if the trait type in the trait table is correct
396 if(_Traits[t] != "quant")
397 return error("the selection model \"gaussian\" does not match with the trait type \"%s\"\n", _Traits[t].c_str());
398
399 if( get_parameter("selection_trait_dimension")->isSet() )
400 _selectTraitDimension = (unsigned int)get_parameter_value("selection_trait_dimension");
401 else
403
404 //check if the number of quanti traits modeled match with the dimensionality:
405 if(_popPtr->getTraitPrototype("quant")->get_parameter_value("quanti_traits") != _selectTraitDimension)
406 return error("The number of quantitative traits in \"quanti_traits\" does not match with\
407 the selection trait dimensionality in \"selection_trait_dimension\"\n");
408
409 // set the selection matrices for multi-stage or spatially variable cases
410 if( !setSelectionMatrix() ) return false;
411
412 if(_selectTraitDimension > 1) {
413
414 if(_eVariance > 0)
416 else
418
419 } else {
420
421 if(_eVariance > 0)
423 else
425
426 }
427
428 } else return error("wrong selection model, must be either \"fix\", \"direct\", \"quadratic\", or \"gaussian\".\n");
429
430 }
431
432 if(sel_model != "fix" && !_paramSet->isSet("selection_trait"))
433 return error("trait under selection is not set, please add parameter \"selection_trait\"\n");
434
435 return true;
436}
TraitPrototype * getTraitPrototype(trait_t type)
Accessor to a TraitPrototype.
Definition: indfactory.cc:140
double getFitnessUnivariateGaussian_VE(Individual *offsprg, unsigned int patch, unsigned int trait)
Returns the fitness of an individual following the Gaussian selection model with several traits under...
Definition: LCEselection.cc:742
bool setSelectionMatrix()
Definition: LCEselection.cc:440
double getFitnessFixedEffect(Individual *offsprg, unsigned int patch, unsigned int trait)
Returns the fitness of an individual in the fixed selection model.
Definition: LCEselection.h:196
double getFitnessDirect(Individual *offsprg, unsigned int patch, unsigned int trait)
Returns the fitness of an individual following the direct selection model.
Definition: LCEselection.h:202
double getFitnessUnivariateQuadratic(Individual *ind, unsigned int patch, unsigned int trait)
Quadratic fitness surface, approximates the Gaussian model for weak selection and/or small deviation ...
Definition: LCEselection.cc:690
double getFitnessMultivariateGaussian_VE(Individual *offsprg, unsigned int patch, unsigned int trait)
Returns the fitness of an individual following the Gaussian selection model with one trait under sele...
Definition: LCEselection.cc:758
double getFitnessMultivariateGaussian(Individual *offsprg, unsigned int patch, unsigned int trait)
Returns the fitness of an individual following the Gaussian selection model with one trait under sele...
Definition: LCEselection.cc:720
double _Fpedigree[5]
Array of pedigree values used in the fixed-selection model.
Definition: LCEselection.h:72
double _letheq
Definition: LCEselection.h:88
vector< string > _SelectionModels
The selection models associated with each trait under selection.
Definition: LCEselection.h:114
vector< string > _Traits
The list of trait types under selection.
Definition: LCEselection.h:108
double getFitnessUnivariateGaussian(Individual *offsprg, unsigned int patch, unsigned int trait)
Returns the fitness of an individual following the Gaussian selection model with several traits under...
Definition: LCEselection.cc:705
vector< string > getMultiArgs()
Definition: param.cc:122
unsigned int getNbCols()
Gives the number of columns.
Definition: tmatrix.h:169
void warning(const char *str,...)
Definition: output.cc:58

References _base_fitness, _eVariance, _FitnessFixModel, _Fpedigree, _getRawFitness, _letheq, SimComponent::_paramSet, LifeCycleEvent::_popPtr, _SelectionModels, _selectTraitDimension, _TraitIndices, _Traits, error(), TMatrix::get(), SimComponent::get_parameter(), SimComponent::get_parameter_value(), getFitnessDirect(), getFitnessFixedEffect(), getFitnessMultivariateGaussian(), getFitnessMultivariateGaussian_VE(), getFitnessUnivariateGaussian(), getFitnessUnivariateGaussian_VE(), getFitnessUnivariateQuadratic(), Param::getMatrix(), Param::getMultiArgs(), TMatrix::getNbCols(), IndFactory::getTraitPrototype(), ParamSet::isSet(), setSelectionMatrix(), and warning().

Referenced by LCE_Selection_base(), and setParameters().

+ Here is the caller graph for this function:

◆ set_std_rate_of_change()

void LCE_Selection_base::set_std_rate_of_change ( )
1000{
1002
1003 //reset rate_of_change matrix in each new replicate:
1004 TMatrix tmpMat;
1005
1006 if(!get_parameter("selection_std_rate_environmental_change")->isMatrix()) {
1007
1008 double val = get_parameter_value("selection_std_rate_environmental_change");
1009
1011 tmpMat.assign(val);
1012
1013 } else {
1014 get_parameter("selection_std_rate_environmental_change")->getMatrix(&tmpMat);
1015 }
1016
1017 unsigned int nCol = tmpMat.ncols(), nRow = tmpMat.nrows();
1018
1020
1021 //copy values, with pattern propagation
1022 for (int p = 0; p < _popPtr->getPatchNbr(); ++p) {
1023 for (int i = 0; i < _selectTraitDimension; ++i) {
1024 _rate_of_change_local_optima.set(p, i, tmpMat.get(p % nRow, i % nCol));
1025 }
1026 }
1027
1028 double* SD = new double [_selectTraitDimension];
1029
1030 for (int i = 0; i < _selectTraitDimension; ++i) {
1031 SD[i] = 0;
1032 }
1033
1034 // check if SD is taken in a reference patch instead of the whole pop
1035 if (_std_rate_reference_patch > -1) {
1036
1039
1040 } else {
1041 // compute SD as the mean within-patch SD
1042
1043 unsigned int cnt = 0;
1044
1045 //get SD only in extant demes, to avoid nans
1046 for(unsigned int i = 0; i < _popPtr->getPatchNbr(); ++i) {
1047 if (_popPtr->getPatch(i)->size(OFFSx) != 0 ) {
1048 cnt++;
1049 addPhenotypicSD(i, SD);
1050 }
1051 }
1052
1053 //compute mean within-patch phenotypic standard deviation:
1054 for (int i = 0; i < _selectTraitDimension; ++i) {
1055 SD[i]/= cnt;
1056
1057 }
1058 } //end if
1059
1060 //multiply the per-trait rates of change by SD:
1061 for (int p = 0; p < _popPtr->getPatchNbr(); ++p) {
1062 for (int i = 0; i < _selectTraitDimension; ++i){
1064 }
1065
1066 }
1067 //log the rates in the simulation log file:
1068 SIMenv::MainSim->_FileServices.log("#selection_rate_environmental_change " +
1070
1071 //compute the change of local optima for current generation:
1073
1074 //now reset the function pointer to changeLocalOptima() for next generation:
1076
1077 delete [] SD;
1078 }
1079}
void log(string message)
Write to the parameter logfile.
Definition: fileservices.cc:359
void addPhenotypicSD(unsigned int deme, double *stDev)
Definition: LCEselection.cc:1083
static SimRunner * MainSim
Definition: simenv.h:42
FileServices _FileServices
Definition: simulation.h:100
void multi(unsigned int i, unsigned int j, double value)
Multiply an element of the matrix by a value.
Definition: tmatrix.h:235
string to_string()
Writes the matrix into a string in Nemo's matrix input format.
Definition: tmatrix.h:292

References SimRunner::_FileServices, LifeCycleEvent::_popPtr, _rate_of_change_local_optima, _selectTraitDimension, _set_std_rate_at_generation, _setNewLocalOptima, _std_rate_reference_patch, addPhenotypicSD(), TMatrix::assign(), changeLocalOptima(), TMatrix::get(), SimComponent::get_parameter(), SimComponent::get_parameter_value(), Metapop::getCurrentGeneration(), Param::getMatrix(), Metapop::getPatch(), Metapop::getPatchNbr(), FileServices::log(), SIMenv::MainSim, TMatrix::multi(), TMatrix::ncols(), TMatrix::nrows(), OFFSx, TMatrix::reset(), TMatrix::set(), Patch::size(), and TMatrix::to_string().

Referenced by execute(), and set_param_rate_of_change().

+ Here is the caller graph for this function:

◆ setLocalOptima()

bool LCE_Selection_base::setLocalOptima ( )

◆ setMeanFitness()

double LCE_Selection_base::setMeanFitness ( age_idx  age)
inline

Sets the _mean_fitness variable to the value of the mean population fitness.

Parameters
agethe age-class index
184{ return (_mean_fitness = getMeanFitness(age));}
double getMeanFitness(age_idx age)
Computes the mean fitness of the whole population for a given age class.
Definition: LCEselection.cc:791

References _mean_fitness, and getMeanFitness().

◆ setMeans()

void LCE_Selection_base::setMeans ( unsigned int  tot_ind)

Computes the average fitness of each pedigree class.

Called after selection has been done in the whole population.

Parameters
tot_indcount of the total number of individuals produced
984{
985 _mean_fitness = 0;
986
987 for(unsigned int i = 0; i < 5; i++) {
989 _fitness[i] /= _ind_cntr[i];
990 _survival[i] /= _ind_cntr[i];
991 _ind_cntr[i] /= tot_ind;
992 }
993
994 _mean_fitness /= tot_ind;
995}

References _fitness, _ind_cntr, _mean_fitness, and _survival.

Referenced by LCE_Breed_Selection::execute(), LCE_Breed_Selection_Disperse::execute(), and execute().

+ Here is the caller graph for this function:

◆ setParameters()

bool LCE_Selection_base::setParameters ( )
virtual

Implements SimComponent.

Reimplemented in LCE_Breed_Selection, and LCE_Breed_Selection_Disperse.

153{
154 //selection may act on more than one trait at a time
155 _Traits.clear();
156 _Traits = get_parameter("selection_trait")->getMultiArgs();
157
158 _TraitIndices.clear();
159
160 for(unsigned int i = 0; i < _Traits.size(); i++) {
161
162 if(_popPtr->getTraitIndex(_Traits[i].c_str()) == -1) {
163
164 return error("cannot attach trait \"%s\" to life cycle event \"%s\", trait has not been initiated.\n",
165 _Traits[i].c_str(), _event_name.c_str());
166
167 } else {
168 _TraitIndices.push_back(_popPtr->getTraitIndex(_Traits[i].c_str()));
169
170
171#ifdef _DEBUG_
172 cout<<"#### attaching trait \""<<_Traits[i].c_str()<<"\" to selection LCE, index = "<<_TraitIndices[i]<<endl;
173#endif
174
175 }
176
177 }
178
179
180 if(!set_fit_model()) return false;
181
182 if(!set_sel_model()) return false;
183
184 if(!set_local_optima()) return false;
185
186 if(!set_param_rate_of_change()) return false;
187
189 _scaling_factor = 1;
191
192 assert(_TraitIndices.size() == _getRawFitness.size());
193
194#ifdef _DEBUG_
195 cout<<"#### selection is acting on "<<_TraitIndices.size()<<" traits\n";
196 cout<<"#### selection is using "<< _getRawFitness.size() <<" fitness functions\n";
197#endif
198
199
200 return true;
201}
int getTraitIndex(trait_t type)
Gives the index of trait with type.
Definition: indfactory.cc:128
std::string _event_name
The param name to be read in the init file.
Definition: lifecycleevent.h:77

References LifeCycleEvent::_event_name, _getRawFitness, _max_fitness, _mean_fitness, LifeCycleEvent::_popPtr, _scaling_factor, _TraitIndices, _Traits, error(), SimComponent::get_parameter(), Param::getMultiArgs(), IndFactory::getTraitIndex(), resetCounters(), set_fit_model(), set_local_optima(), set_param_rate_of_change(), and set_sel_model().

Referenced by LCE_Breed_Selection::setParameters(), and LCE_Breed_Selection_Disperse::setParameters().

+ Here is the caller graph for this function:

◆ setScalingFactorAbsolute()

void LCE_Selection_base::setScalingFactorAbsolute ( age_idx  age,
unsigned int  p 
)
inline

Resets the fitness scaling factor equal to one.

Parameters
agethe age-class index
pthe focal patch index
267{ _scaling_factor = 1; }

References _scaling_factor.

Referenced by set_fit_model().

+ Here is the caller graph for this function:

◆ setScalingFactorGlobal()

void LCE_Selection_base::setScalingFactorGlobal ( age_idx  age,
unsigned int  p 
)

Sets the fitness scaling factor equal to the inverse of the mean population fitness.

Function exits on p != 0; the mean population fitness is computed only once in the execute() procedure, at the beginning of the patch loop.

Parameters
agethe age-class index
pthe focal patch index
870{
871 if(p != 0) return; //we compupte the mean fitness only once
872 _scaling_factor = 1;
874}

References _scaling_factor, and getMeanFitness().

Referenced by set_fit_model(), and LCE_SelectionSH::setDataTable().

+ Here is the caller graph for this function:

◆ setScalingFactorLocal()

void LCE_Selection_base::setScalingFactorLocal ( age_idx  age,
unsigned int  p 
)

Sets the fitness scaling factor equal to the inverse of the mean local patch fitness.

Parameters
agethe age-class index
pthe focal patch index
862{
863 _scaling_factor = 1; //this to have the raw mean fitness below
865}
double getMeanPatchFitness(age_idx age, unsigned int p)
Computes the mean fitness in a given patch for a given age class.
Definition: LCEselection.cc:809

References _scaling_factor, and getMeanPatchFitness().

Referenced by set_fit_model(), and LCE_SelectionSH::setDataTable().

+ Here is the caller graph for this function:

◆ setScalingFactorMaxGlobal()

void LCE_Selection_base::setScalingFactorMaxGlobal ( age_idx  age,
unsigned int  p 
)

Sets the fitness scaling factor equal to the inverse of the maximum population fitness value.

Function exits on p != 0; the mean population fitness is computed only once in the execute() procedure, at the beginning of the patch loop.

Parameters
agethe age-class index
pthe focal patch index
888{
889 if(p != 0) return; //we compute the max fitness only once
890 _scaling_factor = 1;
892}
double getMaxFitness(age_idx age)
Computes the maximum fitness value of the whole population for a given age class.
Definition: LCEselection.cc:825

References _scaling_factor, and getMaxFitness().

Referenced by set_fit_model().

+ Here is the caller graph for this function:

◆ setScalingFactorMaxLocal()

void LCE_Selection_base::setScalingFactorMaxLocal ( age_idx  age,
unsigned int  p 
)

Sets the fitness scaling factor equal to the inverse of the maximum local patch fitness value.

Parameters
agethe age-class index
pthe focal patch index
880{
881 _scaling_factor = 1; //this to have the raw fitness when computing fitnesses below
883}

References _scaling_factor, and getMaxPatchFitness().

Referenced by set_fit_model().

+ Here is the caller graph for this function:

◆ setSelectionMatrix()

bool LCE_Selection_base::setSelectionMatrix ( )
441{
442 TMatrix tmp_mat;
443 unsigned int patchNbr = _popPtr->getPatchNbr();
444
445 if(!get_parameter("selection_matrix")->isSet() && !get_parameter("selection_variance")->isSet())
446 return error("\"selection_matrix\" or \"selection_variance\" must be set with selection model = \"gaussian\".\n");
447
448
449 if(get_parameter("selection_variance")->isSet() && !get_parameter("selection_trait_dimension")->isSet())
450 return error("parameter \"selection_trait_dimension\" is missing!\n");
451
452
453 //clear the selection matrix container
454 vector< TMatrix* >::iterator selIT = _selection_matrix.begin();
455
456 for(; selIT != _selection_matrix.end(); ++selIT) if((*selIT)) delete (*selIT);
457
458 _selection_matrix.clear();
459
460
461 if(get_parameter("selection_matrix")->isSet()) {
462
463 //selection matrix provided, same selection surface in each patch
464 _paramSet->getMatrix("selection_matrix", &tmp_mat);
465
466 if(tmp_mat.getNbCols() != tmp_mat.getNbRows())
467 return error("\"selection_matrix\" must be a square matrix!\n");
468
470
471 //we have one selection matrix per patch, copy it in the container for each patch
472 for(unsigned int i = 0; i < patchNbr; ++i)
473 _selection_matrix.push_back( new TMatrix(tmp_mat) );
474
475 } else {
476 //we have to check for spatial variation in variance and covariances
477 _selectTraitDimension = (unsigned int)get_parameter_value("selection_trait_dimension");
478
479 TMatrix var_spatmat, corr_spatmat;
480
481 //setting variance spatial matrix:
482 var_spatmat.reset(patchNbr, (unsigned)_selectTraitDimension);
483
484 if(get_parameter("selection_variance")->isMatrix()) {
485
486 _paramSet->getMatrix("selection_variance", &tmp_mat);
487
488 if( !setSpatialMatrix("selection_variance","\"selection_trait_dimension\"", &tmp_mat, &var_spatmat,
489 (unsigned)_selectTraitDimension, patchNbr, _paramSet->isSet("selection_randomize") ) )
490 return false;
491
492 } else {
493
494 var_spatmat.assign(get_parameter_value("selection_variance"));
495 }
496
497 //setting correlation spatial matrix:
498 corr_spatmat.reset(patchNbr, (unsigned)_selectTraitDimension*(_selectTraitDimension-1)/2);
499
500 if(get_parameter("selection_correlation")->isMatrix()) {
501
502 _paramSet->getMatrix("selection_correlation", &tmp_mat);
503
504 if( !setSpatialMatrix("selection_correlation","the num of correlation coefficients", &tmp_mat, &corr_spatmat,
505 (unsigned)_selectTraitDimension*(_selectTraitDimension-1)/2, patchNbr, _paramSet->isSet("selection_randomize") ) )
506 return false;
507
508 } else {
509
510 corr_spatmat.assign((get_parameter("selection_correlation")->isSet() ?
511 get_parameter_value("selection_correlation") : 0.0 ));
512 }
513
514 //set the selection matrix:
516 double covar;
517 unsigned int col;
518 for( unsigned int p = 0; p < patchNbr; p++) {
519 col = 0;
520 for( int i = 0; i < _selectTraitDimension; i++) {
521 tmp_mat.set(i, i, var_spatmat.get(p, i));
522 for( int j = i+1; j < _selectTraitDimension; j++) {
523 covar = corr_spatmat.get(p, col) * sqrt( var_spatmat.get(p, i) * var_spatmat.get(p, j) );
524 tmp_mat.set(i, j, covar);
525 tmp_mat.set(j, i, covar);
526 col++;
527 }
528 }
529 _selection_matrix.push_back(new TMatrix(tmp_mat));
530 }
531 }
532
533 if(_selectTraitDimension > 1) {
534
535 //selection on more than one trait
536
537 //inverting the selection matrices:
538 if(_gsl_selection_matrix.size() != 0)
539 for(unsigned int i = 0; i < _gsl_selection_matrix.size(); ++i)
540 if(_gsl_selection_matrix[i] != NULL) gsl_matrix_free( _gsl_selection_matrix[i] );
541
542 _gsl_selection_matrix.clear();
543
544 for( unsigned int p = 0; p < patchNbr; p++) {
545 _selection_matrix[p]->inverse();
546
548
549 _selection_matrix[p]->get_gsl_matrix(_gsl_selection_matrix[p]);
550 }
551
552 //allocate the vectors used by the fitness function:
553 if(_diffs != NULL) gsl_vector_free(_diffs);
554 _diffs = gsl_vector_alloc( _selectTraitDimension );
555
556 if(_res1 != NULL) gsl_vector_free(_res1);
557 _res1 = gsl_vector_alloc( _selectTraitDimension );
558
559 } else {
560
561 //selection on one trait only
562 _selection_variance.clear();
563
564 for(unsigned int p = 0; p < patchNbr; p++)
565 _selection_variance.push_back( _selection_matrix[p]->get(0, 0) );
566 }
567
568
569 return true;
570}
unsigned int getNbRows()
Gives the number of rows.
Definition: tmatrix.h:166

References _diffs, _gsl_selection_matrix, SimComponent::_paramSet, LifeCycleEvent::_popPtr, _res1, _selection_matrix, _selection_variance, _selectTraitDimension, TMatrix::assign(), error(), TMatrix::get(), SimComponent::get_parameter(), SimComponent::get_parameter_value(), ParamSet::getMatrix(), TMatrix::getNbCols(), TMatrix::getNbRows(), Metapop::getPatchNbr(), Param::isMatrix(), Param::isSet(), ParamSet::isSet(), TMatrix::reset(), TMatrix::set(), and setSpatialMatrix().

Referenced by set_sel_model().

+ Here is the caller graph for this function:

Friends And Related Function Documentation

◆ LCE_Breed_Selection

friend class LCE_Breed_Selection
friend

◆ LCE_SelectionFH

friend class LCE_SelectionFH
friend

Referenced by loadFileServices().

◆ LCE_SelectionSH

friend class LCE_SelectionSH
friend

The StatHandler associated class is a friend.

Referenced by loadStatServices().

Member Data Documentation

◆ _base_fitness

◆ _diffs

gsl_vector* LCE_Selection_base::_diffs
private

◆ _do_change_local_opt

bool LCE_Selection_base::_do_change_local_opt
private

◆ _eVariance

double LCE_Selection_base::_eVariance
protected

◆ _fitness

double LCE_Selection_base::_fitness[5]
protected

◆ _FitnessFixModel

double LCE_Selection_base::_FitnessFixModel[5]
private

Absolute fitness values of the five pedigree class for the fixed selection model (lethal equivalents model).

Fitness is computed as: _FitnessFixModel[i] = _base_fitness * exp( -_letheq * _Fpedigree[i] );

Referenced by getFitnessFixedEffect(), and set_sel_model().

◆ _Fpedigree

double LCE_Selection_base::_Fpedigree[5]
private

Array of pedigree values used in the fixed-selection model.

Referenced by set_sel_model().

◆ _getFitness

double(LCE_Selection_base::* LCE_Selection_base::_getFitness) (Individual *, unsigned int)
protected

Pointer to the function returning the individual fitness.

May point to LCE_Selection_base::getFitnessAbsolute or LCE_Selection_base::getFitnessRelative.

Referenced by getFitness(), and set_fit_model().

◆ _getRawFitness

vector< double (LCE_Selection_base::* ) (Individual*, unsigned int, unsigned int) > LCE_Selection_base::_getRawFitness
protected

A vector containing pointers to fitness function related to each trait under selection.

Referenced by getFitnessAbsolute(), LCE_SelectionFH::print(), set_sel_model(), and setParameters().

◆ _gsl_selection_matrix

vector< gsl_matrix* > LCE_Selection_base::_gsl_selection_matrix
private

◆ _ind_cntr

◆ _is_absolute

bool LCE_Selection_base::_is_absolute
protected

◆ _is_local

bool LCE_Selection_base::_is_local
protected

◆ _letheq

double LCE_Selection_base::_letheq
protected

Referenced by set_sel_model().

◆ _local_optima

◆ _max_fitness

double LCE_Selection_base::_max_fitness
protected

Referenced by setParameters().

◆ _mean_fitness

◆ _phe

double* LCE_Selection_base::_phe
private

◆ _rate_of_change_is_std

bool LCE_Selection_base::_rate_of_change_is_std
private

◆ _rate_of_change_local_optima

TMatrix LCE_Selection_base::_rate_of_change_local_optima
private

◆ _res1

gsl_vector * LCE_Selection_base::_res1
private

◆ _scaling_factor

◆ _selection_matrix

vector< TMatrix* > LCE_Selection_base::_selection_matrix
private

◆ _selection_variance

vector< double > LCE_Selection_base::_selection_variance
protected

◆ _SelectionModels

vector< string > LCE_Selection_base::_SelectionModels
protected

The selection models associated with each trait under selection.

Referenced by set_sel_model().

◆ _selectTraitDimension

◆ _set_std_rate_at_generation

unsigned int LCE_Selection_base::_set_std_rate_at_generation
private

◆ _setNewLocalOptima

void(LCE_Selection_base::* LCE_Selection_base::_setNewLocalOptima) (void)
protected

Pointer to the function used to change the local phenotypic optima or its rate of change.

May point to LCE_Selection_base::changeLocalOptima or LCE_Selection_base::set_std_rate_of_change.

Referenced by execute(), set_param_rate_of_change(), and set_std_rate_of_change().

◆ _setScalingFactor

void(LCE_Selection_base::* LCE_Selection_base::_setScalingFactor) (age_idx, unsigned int)
protected

Pointer to the function used to set the fitness scaling factor when fitness is relative.

May point to LCE_Selection_base::setScalingFactorGlobal, LCE_Selection_base::setScalingFactorLocal, or LCE_Selection_base::setScalingFactorAbsolute (returns 1).

Referenced by LCE_Breed_Selection::do_breed_selection_FecFitness(), execute(), and set_fit_model().

◆ _stater

LCE_SelectionSH* LCE_Selection_base::_stater
protected

◆ _std_rate_reference_patch

int LCE_Selection_base::_std_rate_reference_patch
private

◆ _survival

◆ _TraitIndices

vector< unsigned int > LCE_Selection_base::_TraitIndices
protected

◆ _Traits

vector< string > LCE_Selection_base::_Traits
protected

The list of trait types under selection.

Referenced by set_sel_model(), and setParameters().

◆ _writer

LCE_SelectionFH* LCE_Selection_base::_writer
protected

Referenced by loadFileServices().


The documentation for this class was generated from the following files:

Generated for Nemo v2.3.56 by  doxygen 1.9.0 -- Nemo is hosted on  Download Nemo

Locations of visitors to this page
Catalogued on GSR