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

Prototype class of the bitstring-deleterious mutations trait class. More...

#include <ttdeletmutations_bitstring.h>

+ Inheritance diagram for TProtoDeletMutations_bitstring:
+ Collaboration diagram for TProtoDeletMutations_bitstring:

Public Member Functions

 TProtoDeletMutations_bitstring ()
 
 TProtoDeletMutations_bitstring (const TProtoDeletMutations_bitstring &T)
 
 ~TProtoDeletMutations_bitstring ()
 
Getters
int get_nb_locus ()
 
double get_mut_rate ()
 
double get_backmutation_rate ()
 
double get_strength ()
 
double get_dominance ()
 
int get_dominance_model ()
 
bool get_iscontinuous ()
 
double get_init_freq ()
 
double get_fitness_scaling_factor ()
 
float * get_s_continous ()
 
float * get_hs_continous ()
 
float ** get_effects () const
 
float get_effect (unsigned int at, unsigned int loc)
 
Setters
void set_effects ()
 
void set_effect (unsigned int at, unsigned int loc, float val)
 
double set_effects_exp ()
 
double set_effects_gamma ()
 
double set_effects_lognorm ()
 
void reset_effect_table ()
 
void write_effects_to_parameter ()
 
Parameter setters
bool setSelectionParameters ()
 
bool setEffectsFromInput ()
 
Inheritance functions
void inherit_low (sex_t SEX, bitstring *seq, bitstring **parent)
 
void inherit_free (sex_t SEX, bitstring *seq, bitstring **parent)
 
TraitPrototype implementations
virtual void init ()
 
virtual void reset ()
 
virtual TTDeletMutations_bitstringhatch ()
 
virtual TProtoDeletMutations_bitstringclone ()
 
virtual trait_t get_type () const
 
StorageComponent implementation
virtual void store_data (BinaryStorageBuffer *saver)
 
virtual bool retrieve_data (BinaryStorageBuffer *reader)
 
SimComponent implementation
virtual bool setParameters ()
 
virtual void loadFileServices (FileServices *loader)
 
virtual void loadStatServices (StatServices *loader)
 
virtual bool resetParameterFromSource (std::string param, SimComponent *cmpt)
 
- Public Member Functions inherited from TTProtoWithMap
 TTProtoWithMap ()
 
 TTProtoWithMap (const TTProtoWithMap &TP)
 
virtual ~TTProtoWithMap ()
 
void setMapIndex (unsigned int idx)
 
unsigned int getMapIndex ()
 
bool setGeneticMapParameters (string prefix)
 
void addGeneticMapParameters (string prefix)
 
bool setRecombinationMapRandom ()
 
bool setRecombinationMapNonRandom (vector< vector< double > > *lociPositions)
 
bool setRecombinationMapFixed ()
 
bool setNumLociPerChromosome (string param_name)
 
void reset_recombination_pointers ()
 
void registerGeneticMap ()
 
void unregisterFromGeneticMap ()
 
bool areGeneticMapParamSet (string prefix)
 
bool isRecombinationFree (string prefix)
 
void recordRandomMap ()
 
virtual void reset ()
 
- Public Member Functions inherited from TraitPrototype
virtual void reset ()=0
 Called at the end of a simulation to reset the Traits' prototypes (e.g. More...
 
virtual TTraithatch ()=0
 Creates the trait of which it is the prototype, called by IndFactory::makePrototype(). More...
 
virtual TraitPrototypeclone ()=0
 Returns a copy of itself. More...
 
virtual trait_t get_type () const =0
 Type accessor. More...
 
virtual void set_index (int idx)
 Sets the traits index. More...
 
virtual int get_index ()
 Index getter. More...
 
- Public Member Functions inherited from StorableComponent
virtual void store_data (BinaryStorageBuffer *saver)=0
 Interface to store the component data (e.g. gene values) into a binary buffer. More...
 
virtual bool retrieve_data (BinaryStorageBuffer *reader)=0
 Interface to retrieve the same data from the binary buffer. More...
 
virtual ~StorableComponent ()
 
- 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_paramset (std::string name, bool required, SimComponent *owner)
 Sets a new ParamSet and name it. 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...
 

Private Attributes

unsigned int _nb_locus
 
unsigned int _fitness_model
 
unsigned int _mutation_model
 
int _dominance_model
 
double _fitness_scaling_factor
 
double _init_freq
 
double _mut_rate
 
double _back_mutation_rate
 
double _strength
 
double _dominance
 
double _dist_p1
 
double _dist_p2
 
bool _continuous_effects
 
double(TTDeletMutations_bitstring::* _viability_func_ptr )(void)
 
double(TProtoDeletMutations_bitstring::* _set_effects_func )(void)
 
void(TProtoDeletMutations_bitstring::* _inherit_func_ptr )(sex_t, bitstring *, bitstring **)
 
TTDeletMutBitstrSH_stats
 
TTDeletMutBitstrFH_writer
 
TTDeletMutBitstrFH_reader
 
float ** _effects
 

Additional Inherited Members

- Static Public Member Functions inherited from TTProtoWithMap
static void recombine (unsigned long indID)
 
- Static Public Attributes inherited from TTProtoWithMap
static GeneticMap _map
 
- Protected Attributes inherited from TTProtoWithMap
unsigned int _mapIndex
 
double _totRecombEventsMean
 
double _recombRate
 
double _mapResolution
 
unsigned int _numChromosome
 
unsigned int _numLoci
 
double * _recombRatePerChrmsm
 
unsigned int_numLociPerChrmsm
 
unsigned int_chrsmLength
 
unsigned int_lociMapPositions
 
- Protected Attributes inherited from TraitPrototype
int _index
 The trait index in the Individual traits table. More...
 
- Protected Attributes inherited from SimComponent
ParamSet_paramSet
 The parameters container. More...
 

Detailed Description

Prototype class of the bitstring-deleterious mutations trait class.

Constructor & Destructor Documentation

◆ TProtoDeletMutations_bitstring() [1/2]

TProtoDeletMutations_bitstring::TProtoDeletMutations_bitstring ( )
54_writer(0), _reader(0), _effects(0)
55{
56 set_paramset("delet", false, this);
57
58 add_parameter("delet_loci",INT,true,false,0,0);
59
60 add_parameter("delet_mutation_rate",DBL,true,true,0,1, 0);
61 add_parameter("delet_backmutation_rate",DBL,false,true,0,1, 0);
62 add_parameter("delet_mutation_model",INT,false,true,1,2, 0);
63 add_parameter("delet_init_freq",DBL,false,true,0,1, 0);
64
65 //genetic map parameters:
67
68 //effects parameters
69 add_parameter("delet_effects_mean",DBL,false,false,0,1, 0);
70 add_parameter("delet_dominance_mean",DBL,false,true,0,1, 0);
71 add_parameter("delet_fitness_model",INT,false,true,1,2, 0);
72 add_parameter("delet_effects", MAT, false, false, 0, 0, 0);
73 add_parameter("delet_effects_distribution",STR,false,false,0,0, 0);
74 add_parameter("delet_effects_dist_param1",DBL,false,false,0,0, 0);
75 add_parameter("delet_effects_dist_param2",DBL,false,false,0,0, 0);
76 add_parameter("delet_fitness_scaling_factor",DBL,false,false,0,0, 0);
77 add_parameter("delet_dominance_model", INT, false, true, 1, 2, 0);
78 //for back compatibility:
79// add_parameter("delet_sel_coef",DBL,false,true,0,1, 0);
80// add_parameter("delet_dom_coef",DBL,false,true,0,1, 0);
81 //output parameters:
82 add_parameter("delet_save_genotype",BOOL,false,false,0,0);
83 add_parameter("delet_genot_dir",STR,false,false,0,0);
84 add_parameter("delet_genot_logtime",INT,false,false,0,0);
85
86}
virtual void set_paramset(ParamSet *paramset)
Sets the ParamSet member.
Definition: simcomponent.h:86
virtual void add_parameter(Param *param)
Interface to add a parameter to the set.
Definition: simcomponent.h:112
TTDeletMutBitstrFH * _writer
Definition: ttdeletmutations_bitstring.h:266
unsigned int _mutation_model
Definition: ttdeletmutations_bitstring.h:250
double _strength
Definition: ttdeletmutations_bitstring.h:256
unsigned int _nb_locus
Definition: ttdeletmutations_bitstring.h:248
float ** _effects
Definition: ttdeletmutations_bitstring.h:268
double _fitness_scaling_factor
Definition: ttdeletmutations_bitstring.h:252
void(TProtoDeletMutations_bitstring::* _inherit_func_ptr)(sex_t, bitstring *, bitstring **)
Definition: ttdeletmutations_bitstring.h:263
double _init_freq
Definition: ttdeletmutations_bitstring.h:253
double _mut_rate
Definition: ttdeletmutations_bitstring.h:254
double _dominance
Definition: ttdeletmutations_bitstring.h:257
int _dominance_model
Definition: ttdeletmutations_bitstring.h:251
double _back_mutation_rate
Definition: ttdeletmutations_bitstring.h:255
unsigned int _fitness_model
Definition: ttdeletmutations_bitstring.h:249
TTDeletMutBitstrFH * _reader
Definition: ttdeletmutations_bitstring.h:267
double(TTDeletMutations_bitstring::* _viability_func_ptr)(void)
Definition: ttdeletmutations_bitstring.h:261
bool _continuous_effects
Definition: ttdeletmutations_bitstring.h:260
TTDeletMutBitstrSH * _stats
Definition: ttdeletmutations_bitstring.h:265
void addGeneticMapParameters(string prefix)
Definition: ttrait_with_map.cc:77
@ 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(), TTProtoWithMap::addGeneticMapParameters(), BOOL, DBL, INT, MAT, SimComponent::set_paramset(), and STR.

Referenced by clone().

+ Here is the caller graph for this function:

◆ TProtoDeletMutations_bitstring() [2/2]

TProtoDeletMutations_bitstring::TProtoDeletMutations_bitstring ( const TProtoDeletMutations_bitstring T)

◆ ~TProtoDeletMutations_bitstring()

TProtoDeletMutations_bitstring::~TProtoDeletMutations_bitstring ( )
104{
105 if(_stats != NULL) {delete _stats; _stats = NULL;}
106 if(_reader) {delete _reader; _reader= NULL;}
107 if(_writer != NULL) {delete _writer; _writer = NULL;}
108 if(_effects != NULL) {
109 delete [] _effects[0];
110 delete [] _effects[1];
111 delete [] _effects;
112 _effects = NULL;
113 }
114
115}

References _effects, _reader, _stats, and _writer.

Member Function Documentation

◆ clone()

virtual TProtoDeletMutations_bitstring * TProtoDeletMutations_bitstring::clone ( )
inlinevirtual

Implements TraitPrototype.

225{return new TProtoDeletMutations_bitstring((*this));}
TProtoDeletMutations_bitstring()
Definition: ttdeletmutations_bitstring.cc:50

References TProtoDeletMutations_bitstring().

◆ get_backmutation_rate()

double TProtoDeletMutations_bitstring::get_backmutation_rate ( )
inline
181{return _back_mutation_rate;}

References _back_mutation_rate.

Referenced by TTDeletMutations_bitstring::mutate_noredraw(), and TTDeletMutations_bitstring::set_mutation_func_ptr().

+ Here is the caller graph for this function:

◆ get_dominance()

double TProtoDeletMutations_bitstring::get_dominance ( )
inline

◆ get_dominance_model()

int TProtoDeletMutations_bitstring::get_dominance_model ( )
inline
184{return _dominance_model;}

References _dominance_model.

Referenced by TTDeletMutations_bitstring::set_allele_value().

+ Here is the caller graph for this function:

◆ get_effect()

float TProtoDeletMutations_bitstring::get_effect ( unsigned int  at,
unsigned int  loc 
)
inline

◆ get_effects()

float ** TProtoDeletMutations_bitstring::get_effects ( ) const
inline
190{return _effects;}

References _effects.

Referenced by resetParameterFromSource().

+ Here is the caller graph for this function:

◆ get_fitness_scaling_factor()

double TProtoDeletMutations_bitstring::get_fitness_scaling_factor ( )
inline

◆ get_hs_continous()

float * TProtoDeletMutations_bitstring::get_hs_continous ( )
inline
189{return _effects[0];}

References _effects.

Referenced by TTDeletMutBitstrFH::FHwrite().

+ Here is the caller graph for this function:

◆ get_init_freq()

double TProtoDeletMutations_bitstring::get_init_freq ( )
inline
186{return _init_freq; }

References _init_freq.

Referenced by TTDeletMutations_bitstring::init_sequence().

+ Here is the caller graph for this function:

◆ get_iscontinuous()

bool TProtoDeletMutations_bitstring::get_iscontinuous ( )
inline
185{return _continuous_effects;}

References _continuous_effects.

Referenced by TTDeletMutBitstrFH::FHwrite(), and TTDeletMutBitstrSH::TTDeletMutBitstrSH().

+ Here is the caller graph for this function:

◆ get_mut_rate()

double TProtoDeletMutations_bitstring::get_mut_rate ( )
inline
180{return _mut_rate;}

References _mut_rate.

◆ get_nb_locus()

int TProtoDeletMutations_bitstring::get_nb_locus ( )
inline

◆ get_s_continous()

float * TProtoDeletMutations_bitstring::get_s_continous ( )
inline
188{return _effects[1];}

References _effects.

Referenced by TTDeletMutBitstrFH::FHwrite(), and TTDeletMutBitstrSH::setLethalEquivalents().

+ Here is the caller graph for this function:

◆ get_strength()

double TProtoDeletMutations_bitstring::get_strength ( )
inline

◆ get_type()

virtual trait_t TProtoDeletMutations_bitstring::get_type ( ) const
inlinevirtual

Implements TraitPrototype.

227{return DELE;}
#define DELE
Definition: types.h:66

References DELE.

◆ hatch()

TTDeletMutations_bitstring * TProtoDeletMutations_bitstring::hatch ( )
virtual

Implements TraitPrototype.

385{
387
388 new_trait->set_proto(this);
389
390 new_trait->set_nb_locus(_nb_locus);
391 new_trait->set_mut_rate(_mut_rate,_nb_locus);
392
396
397 return new_trait;
398}
Bitstring implementation of TTDeletMutations with recombination.
Definition: ttdeletmutations_bitstring.h:47
void set_inherit_func_ptr(void(TProtoDeletMutations_bitstring::*theFunc)(sex_t, bitstring *, bitstring **))
Definition: ttdeletmutations_bitstring.h:119
void set_mutation_func_ptr(unsigned int m_model)
Definition: ttdeletmutations_bitstring.cc:617
void set_mut_rate(double val, int nloc)
Definition: ttdeletmutations_bitstring.h:117
void set_viability_func_ptr(unsigned int f_model, bool is_cont)
Definition: ttdeletmutations_bitstring.cc:594
void set_proto(TProtoDeletMutations_bitstring *proto)
Definition: ttdeletmutations_bitstring.h:115
void set_nb_locus(int val)
Definition: ttdeletmutations_bitstring.h:116

References _continuous_effects, _fitness_model, _inherit_func_ptr, _mut_rate, _mutation_model, _nb_locus, TTDeletMutations_bitstring::set_inherit_func_ptr(), TTDeletMutations_bitstring::set_mut_rate(), TTDeletMutations_bitstring::set_mutation_func_ptr(), TTDeletMutations_bitstring::set_nb_locus(), TTDeletMutations_bitstring::set_proto(), and TTDeletMutations_bitstring::set_viability_func_ptr().

◆ inherit_free()

void TProtoDeletMutations_bitstring::inherit_free ( sex_t  SEX,
bitstring seq,
bitstring **  parent 
)
inline
403{
404 for(unsigned int i = 0; i < _nb_locus; )
405 seq[i] = (bool)(*parent[RAND::RandBool()])[i];
406}
static bool RandBool()
Returns a random boolean.
Definition: Uniform.h:163

References _nb_locus, and RAND::RandBool().

◆ inherit_low()

void TProtoDeletMutations_bitstring::inherit_low ( sex_t  SEX,
bitstring seq,
bitstring **  parent 
)
inline
411{
412
413 register unsigned int nbRec, prevLoc = 0, chrm_bloc;
414 register bool flipper;
415
416 vector< unsigned int >& recTable = _map.getRecLoci(SEX, _mapIndex);
417 vector< bool > & firstRecPos = _map.getFirstRecPosition(SEX);
418
419 nbRec = recTable.size();
420
421 // cout << "\nTProtoDeletMutations_bitstring::inherit_low -- "<<nbRec<<" recombination events\n";
422
423 for(unsigned int c = 0, stride = 0, rec = 0; c < _numChromosome; ++c) {
424
425 flipper = firstRecPos[c];
426
427 chrm_bloc = stride + _numLociPerChrmsm[c]; //number of loci considered so far
428
429 prevLoc = stride; //stride is the first locus of a chromosome
430
431 //do recombination chromosome-wise
432 for(; recTable[rec] < chrm_bloc && rec < nbRec; rec++) {
433
434 // cout << " -- copy from "<<prevLoc<<" to "<<recTable[rec]<<" (side="<<flipper<<")"<<endl;
435
436 seq->copy( *parent[flipper], prevLoc, recTable[rec]);
437
438 // memcpy(&seq[prevLoc], &parent[flipper][prevLoc], (recTable[rec] - prevLoc));
439
440 prevLoc = recTable[rec];
441
442 flipper = !flipper;
443 }
444
445 // cout << " -- copy from "<<prevLoc<<" to "<<chrm_bloc<<" (side="<<flipper<<")"<<endl;
446 //copy what's left between the last x-over point and the end of the chrmsme
447 seq->copy( *parent[flipper], prevLoc, chrm_bloc);
448
449 // memcpy(&seq[prevLoc], &parent[flipper][prevLoc], (chrm_bloc - prevLoc));
450
451 stride += _numLociPerChrmsm[c];
452
453 }
454 // cout << "parent chromosomes:\n --0:";
455 // for(unsigned int i = 0; i < _nb_locus; ++i)
456 // cout << (*parent[0])[i];
457 // cout << "\n --1:";
458 // for(unsigned int i = 0; i < _nb_locus; ++i)
459 // cout << (*parent[1])[i];
460 // cout << "\ngamete:\n --0:";
461 // for(unsigned int i = 0; i < _nb_locus; ++i)
462 // cout << (*seq)[i];
463 // cout << "\n";
464}
vector< unsigned int > & getRecLoci(sex_t SEX, unsigned int trait)
Returns a vector of the loci where crossing-overs take place.
Definition: ttrait_with_map.h:157
vector< bool > & getFirstRecPosition(sex_t SEX)
Returns the vector of the first chromosome position for recombination, used for all traits.
Definition: ttrait_with_map.h:163
static GeneticMap _map
Definition: ttrait_with_map.h:201
unsigned int _numChromosome
Definition: ttrait_with_map.h:189
unsigned int * _numLociPerChrmsm
Definition: ttrait_with_map.h:192
unsigned int _mapIndex
Definition: ttrait_with_map.h:185
void copy(const bitstring &b)
Unchecked copy, assumes we have sames sizes.
Definition: bitstring.h:260

References TTProtoWithMap::_map, TTProtoWithMap::_mapIndex, TTProtoWithMap::_numChromosome, TTProtoWithMap::_numLociPerChrmsm, bitstring::copy(), GeneticMap::getFirstRecPosition(), and GeneticMap::getRecLoci().

Referenced by setParameters().

+ Here is the caller graph for this function:

◆ init()

virtual void TProtoDeletMutations_bitstring::init ( )
inlinevirtual
219{setParameters();};
virtual bool setParameters()
Definition: ttdeletmutations_bitstring.cc:119

References setParameters().

◆ loadFileServices()

void TProtoDeletMutations_bitstring::loadFileServices ( FileServices loader)
virtual

Implements SimComponent.

343{
344 // --- THE READER ---
345 //always add the reader:
346 if(_reader) delete _reader;
347 _reader = new TTDeletMutBitstrFH(this);
348 //set to read mode:
350 //attach to file manager:
351 loader->attach_reader(_reader);
352
353 //writer
354 if(get_parameter("delet_save_genotype")->isSet()) {
355
356 if(_writer) delete _writer;
357
358 _writer = new TTDeletMutBitstrFH(this);
359
360 Param* param = get_parameter("delet_genot_logtime");
361
362 if(param->isMatrix()) {
363
364 TMatrix temp;
365 param->getMatrix(&temp);
366 _writer->set_multi(true, true, 1, &temp, get_parameter("delet_genot_dir")->getArg());
367
368 } else {
369 // rpl_per, gen_per, rpl_occ, gen_occ, rank, path, self-ref
370 _writer->set(true, true, 1, (param->isSet() ? (int)param->getValue() : 0),
371 0, get_parameter("delet_genot_dir")->getArg(),this);
372 }
373
374 loader->attach(_writer);
375
376 } else if(_writer != NULL) {
377 delete _writer;
378 _writer = NULL;
379 }
380}
virtual void set_multi(bool rpl_per, bool gen_per, int rpl_occ, TMatrix *Occ, string path)
Definition: filehandler.h:197
void set_isInputHandler(bool val)
Definition: filehandler.h:151
virtual void attach_reader(FileHandler *FH)
Attaches the FileHandler to the current list (_readers) of the FileServices.
Definition: fileservices.cc:73
virtual void attach(Handler *FH)
Attaches the FileHandler to the current list (_writers) of the FileServices.
Definition: fileservices.cc:60
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
A class to handle matrix in params, coerces matrix into a vector of same total size.
Definition: tmatrix.h:49
The FileHandler associated with the TTDeletMutations_bitstring trait.
Definition: ttdeletmutations_bitstring.h:343
virtual void set(bool rpl_per, bool gen_per, int rpl_occ, int gen_occ, int rank, string path, TP *trait_proto)
Definition: filehandler.h:236

References _reader, _writer, FileServices::attach(), FileServices::attach_reader(), SimComponent::get_parameter(), Param::getMatrix(), Param::getValue(), Param::isMatrix(), Param::isSet(), TraitFileHandler< TP >::set(), FileHandler::set_isInputHandler(), and FileHandler::set_multi().

◆ loadStatServices()

void TProtoDeletMutations_bitstring::loadStatServices ( StatServices loader)
virtual

Implements SimComponent.

332{
333 if(_stats != NULL)
334 delete _stats;
335
336 _stats = new TTDeletMutBitstrSH(this);
337
338 loader->attach(_stats);
339}// ----------------------------------------------------------------------------------------
virtual void attach(Handler *H)
attach the StatHandler to the current list (_statHandlers) of the StatServices
Definition: statservices.cc:177
The StatHandler for TTDeletMutations_bitstring.
Definition: ttdeletmutations_bitstring.h:272

References _stats, and StatServices::attach().

◆ reset()

virtual void TProtoDeletMutations_bitstring::reset ( )
inlinevirtual

Reimplemented from TTProtoWithMap.

virtual void reset()
Definition: ttrait_with_map.cc:611

References TTProtoWithMap::reset().

◆ reset_effect_table()

void TProtoDeletMutations_bitstring::reset_effect_table ( )
294{
295 if(_effects != NULL) {
296 delete [] _effects[0];
297 delete [] _effects[1];
298 delete [] _effects;
299 }
300
301 _effects = new float* [2];
302 _effects[0] = new float [_nb_locus];
303 _effects[1] = new float [_nb_locus];
304
305}

References _effects, and _nb_locus.

Referenced by resetParameterFromSource(), set_effects(), and setEffectsFromInput().

+ Here is the caller graph for this function:

◆ resetParameterFromSource()

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

Implements SimComponent.

513{
514 TProtoDeletMutations_bitstring *srceProto = dynamic_cast< TProtoDeletMutations_bitstring* >(cmpt);
515
516 // this works only to reset the fitness effects
517 if(param != "delet_effects") return false;
518
519 float** efx = srceProto->get_effects();
520
521 if(efx == NULL)
522 return error("the binary source population did not properly set the fitness effects of\
523 deleterious mutations, no effects received.\n");
524
525 // the genetic architecture has already been checked in Metapop::loadPopFromBinarySource
526 // we can skip the check of table size
527
528 reset_effect_table(); // deallocate and re-allocate the _effects table
529
530 //copy effects
531 memcpy(_effects[0], efx[0], _nb_locus*sizeof(float));
532 memcpy(_effects[1], efx[1], _nb_locus*sizeof(float));
533
534 //record effects into the parameter
536
537 return true;
538}
Prototype class of the bitstring-deleterious mutations trait class.
Definition: ttdeletmutations_bitstring.h:171
void reset_effect_table()
Definition: ttdeletmutations_bitstring.cc:293
void write_effects_to_parameter()
Definition: ttdeletmutations_bitstring.cc:309
float ** get_effects() const
Definition: ttdeletmutations_bitstring.h:190
int error(const char *str,...)
Definition: output.cc:77

References _effects, _nb_locus, error(), get_effects(), reset_effect_table(), and write_effects_to_parameter().

◆ retrieve_data()

bool TProtoDeletMutations_bitstring::retrieve_data ( BinaryStorageBuffer reader)
virtual

Implements StorableComponent.

483{
484 unsigned int dummy_int;
485 reader->read(&dummy_int,sizeof(int));
486
487 if(dummy_int != _nb_locus ){
488 warning("Deleterious mutation trait::locus number in binary source file differs from parameter value!\n\
489 >>>> adjusting locus number to match source population: %i set to %i\n", _nb_locus, dummy_int);
490 _nb_locus = dummy_int;
491 }
492
493 char dummy_bool = 0;
494 reader->read(&dummy_bool, 1);
495
496 if(dummy_bool != _continuous_effects) {
497 error("Deleterious mutation trait::mutation effects model in binary source file differ from parameter value!\n\
498 >>>> adjusting model to match source population: continuous = %i to continuous = %i\n",
499 _continuous_effects, (int)dummy_bool);
500 _continuous_effects = dummy_bool;
501 }
502
504 reader->read(_effects[0],_nb_locus * sizeof(float));
505 reader->read(_effects[1],_nb_locus * sizeof(float));
506 }
507 return true;
508}
void read(void *out, unsigned int nb_bytes)
Definition: binarystoragebuffer.h:162
void warning(const char *str,...)
Definition: output.cc:58

References _continuous_effects, _effects, _nb_locus, error(), BinaryStorageBuffer::read(), and warning().

◆ set_effect()

void TProtoDeletMutations_bitstring::set_effect ( unsigned int  at,
unsigned int  loc,
float  val 
)
inline
197{_effects[at][loc] = val;}

References _effects.

Referenced by TTDeletMutBitstrFH::FHread(), and TTDeletMutations_bitstring::set_allele_value().

+ Here is the caller graph for this function:

◆ set_effects()

void TProtoDeletMutations_bitstring::set_effects ( )
267{
268 reset_effect_table (); //de-allocate, re-allocate memory
269
270 double dom;
271 double k = -1.0*log(2*_dominance) / _strength;
272
273 for(unsigned int i = 0; i < _nb_locus; ++i){
274
275 do{
276 _effects[1][i] = (float)(this->*_set_effects_func)(); //homozygote effect: s
277 }while(_effects[1][i] > 1); //truncate distribution
278
279 if(_dominance_model == 1)
280 dom = exp(-1.0*_effects[1][i]*k)/2.0; //scaling of h on s
281 else dom = _dominance;
282
283 _effects[0][i] = (float)(dom * _effects[1][i]); //heterozygote effect: hs
284
285 }
286
287 // save the effects to the parameter "delet_effects" will be written to the .log file
289}
double(TProtoDeletMutations_bitstring::* _set_effects_func)(void)
Definition: ttdeletmutations_bitstring.h:262

References _dominance, _dominance_model, _effects, _nb_locus, _set_effects_func, _strength, reset_effect_table(), and write_effects_to_parameter().

Referenced by setSelectionParameters().

+ Here is the caller graph for this function:

◆ set_effects_exp()

double TProtoDeletMutations_bitstring::set_effects_exp ( )
inline
static double Exponential(double mu)
Definition: Uniform.h:448

References _strength, and RAND::Exponential().

Referenced by setSelectionParameters().

+ Here is the caller graph for this function:

◆ set_effects_gamma()

double TProtoDeletMutations_bitstring::set_effects_gamma ( )
inline
199{return RAND::Gamma(_dist_p1, _dist_p2);}
static double Gamma(double a, double b)
Definition: Uniform.h:386
double _dist_p2
Definition: ttdeletmutations_bitstring.h:259
double _dist_p1
Definition: ttdeletmutations_bitstring.h:258

References _dist_p1, _dist_p2, and RAND::Gamma().

Referenced by setSelectionParameters().

+ Here is the caller graph for this function:

◆ set_effects_lognorm()

double TProtoDeletMutations_bitstring::set_effects_lognorm ( )
inline
static double LogNormal(double zeta, double sigma)
Definition: Uniform.h:358

References _dist_p1, _dist_p2, and RAND::LogNormal().

Referenced by setSelectionParameters().

+ Here is the caller graph for this function:

◆ setEffectsFromInput()

bool TProtoDeletMutations_bitstring::setEffectsFromInput ( )
239{
240 TMatrix tmp;
241
242 get_parameter("delet_effects")->getMatrix(&tmp);
243
244 if(tmp.getNbCols() != _nb_locus)
245 return error("the number of mutational effects of the deleterious trait is not equal to its number of loci, please adjust.\n");
246
247 if(tmp.getNbRows() != 2 )
248 return error("the parameter \"delet_effects\" must have two rows, the first for\
249 the mutational effect in heterozygotes, the second for the effects in homozygotes.\n");
250
251 _continuous_effects = true;
252
254
255 for(unsigned int i = 0; i < _nb_locus; ++i){
256
257 _effects[0][i] = (float)tmp.get(0, i); //heterozygote effect: hs
258 _effects[1][i] = (float)tmp.get(1, i); //homozygote effect: s
259 }
260
261 return true;
262}
double get(unsigned int i, unsigned int j)
Accessor to element at row i and column j.
Definition: tmatrix.h:147
unsigned int getNbRows()
Gives the number of rows.
Definition: tmatrix.h:166
unsigned int getNbCols()
Gives the number of columns.
Definition: tmatrix.h:169

References _continuous_effects, _effects, _nb_locus, error(), TMatrix::get(), SimComponent::get_parameter(), Param::getMatrix(), TMatrix::getNbCols(), TMatrix::getNbRows(), and reset_effect_table().

Referenced by setSelectionParameters().

+ Here is the caller graph for this function:

◆ setParameters()

bool TProtoDeletMutations_bitstring::setParameters ( )
virtual

Implements SimComponent.

120{
121 _nb_locus = (int)get_parameter_value("delet_loci");
122 _mut_rate = get_parameter_value("delet_mutation_rate");
123
124
125 //we set the fitness model to 1 (multiplicative by default))
126 if(get_parameter("delet_fitness_model")->isSet())
127 _fitness_model = (int)get_parameter_value("delet_fitness_model");
128 else
129 _fitness_model = 1;
130
131 //reverse mutation rate; we only store the ratio of back/forward mutation rate
132 if( get_parameter("delet_backmutation_rate")->isSet() ) {
133
134 _back_mutation_rate = get_parameter_value("delet_backmutation_rate");
135
137
138 return error("Deleterious mutations: The backward (reverse) mutation rate is larger than\
139the forward mutation rate!\n");
140
141 else
142
144
145
146 } else
147 _back_mutation_rate = 0; //by default we don't process backward mutations
148
149 //mutation model
150 if( get_parameter("delet_mutation_model")->isSet() )
151 _mutation_model = (unsigned int)get_parameter_value("delet_mutation_model");
152 else
153 _mutation_model = 1;
154
155 //initial allele frequency, applies to all loci equally
156 if( (_init_freq = get_parameter_value("delet_init_freq")) == -1.0)
157 _init_freq = 0;
158
159 //selection parameters:
160 if(!setSelectionParameters()) return false;
161
162
163 if( (_fitness_scaling_factor = get_parameter_value("delet_fitness_scaling_factor")) == -1.0)
165
167
169}
virtual double get_parameter_value(std::string name)
Param value getter.
Definition: simcomponent.h:143
void inherit_low(sex_t SEX, bitstring *seq, bitstring **parent)
Definition: ttdeletmutations_bitstring.cc:410
bool setSelectionParameters()
Definition: ttdeletmutations_bitstring.cc:173
bool setGeneticMapParameters(string prefix)
Definition: ttrait_with_map.cc:125

References _back_mutation_rate, _fitness_model, _fitness_scaling_factor, _inherit_func_ptr, _init_freq, _mut_rate, _mutation_model, _nb_locus, error(), SimComponent::get_parameter(), SimComponent::get_parameter_value(), inherit_low(), Param::isSet(), TTProtoWithMap::setGeneticMapParameters(), and setSelectionParameters().

Referenced by init().

+ Here is the caller graph for this function:

◆ setSelectionParameters()

bool TProtoDeletMutations_bitstring::setSelectionParameters ( )
174{
175
176 if(get_parameter("delet_effects")->isSet())
177 return setEffectsFromInput();
178
179 //effects and dominance parameters:
180 _strength = get_parameter_value("delet_effects_mean");
181 _dominance = get_parameter_value("delet_dominance_mean");
182 _dist_p1 = get_parameter_value("delet_effects_dist_param1");
183 _dist_p2 = get_parameter_value("delet_effects_dist_param2");
184
185 //for backward compatibility only, param removed
186 _dominance_model = (int)get_parameter_value("delet_dominance_model");
187
188 if(_dominance_model == -1) _dominance_model = 1; //1 means continuous; 2 means constant
189
190 string distro = _paramSet->getArg("delet_effects_distribution");
191
192 _continuous_effects = true;
193
194 if(distro.empty() || distro == "constant") {
195
196 _continuous_effects = false;
197
199
200 } else if(distro == "exponential") {
201
203
204 } else if(distro == "gamma") {
205
207 //check params:
208 if(_dist_p1 == -1) {
209 return error("missing parameter 1 (shape) of the gamma distribution for the deleterious mutation effects.");
210 }
211
212 if(_dist_p2 == -1)
213 _dist_p2 = _strength / _dist_p1; //mean of gamma == _strength = a*b with a == shape == _dist_p1
214
215 } else if(distro == "lognormal") {
216
218 //check params:
219 if(_dist_p1 == -1) {
220 return error("missing parameter 1 (mu) of the log-normal distribution for the deleterious mutation effects.");
221 }
222 if(_dist_p2 == -1) {
223 return error("missing parameter 2 (sigma) of the log-normal distribution for the deleterious mutation effects.");
224 }
225
226 } else {
227 return error("Deleterious effects distribution \"%s\" not implemented.",distro.c_str());
228 }
229
230
232
233 return true;
234}
string getArg(string name)
Accessor to the parameters argument string.
Definition: param.h:300
double set_effects_lognorm()
Definition: ttdeletmutations_bitstring.h:200
double set_effects_gamma()
Definition: ttdeletmutations_bitstring.h:199
void set_effects()
Definition: ttdeletmutations_bitstring.cc:266
double set_effects_exp()
Definition: ttdeletmutations_bitstring.h:198
bool setEffectsFromInput()
Definition: ttdeletmutations_bitstring.cc:238

References _continuous_effects, _dist_p1, _dist_p2, _dominance, _dominance_model, SimComponent::_paramSet, _set_effects_func, _strength, error(), SimComponent::get_parameter(), SimComponent::get_parameter_value(), ParamSet::getArg(), set_effects(), set_effects_exp(), set_effects_gamma(), set_effects_lognorm(), and setEffectsFromInput().

Referenced by setParameters().

+ Here is the caller graph for this function:

◆ store_data()

void TProtoDeletMutations_bitstring::store_data ( BinaryStorageBuffer saver)
virtual

Implements StorableComponent.

469{
470 saver->store(&_nb_locus,sizeof(int));
471 //store all bool on 1 byte
472 char dummy = _continuous_effects;
473 saver->store(&dummy, 1);
475 saver->store(_effects[0],_nb_locus * sizeof(float));
476 saver->store(_effects[1],_nb_locus * sizeof(float));
477 }
478}
void store(void *stream, unsigned int nb_bytes)
Definition: binarystoragebuffer.cc:16

References _continuous_effects, _effects, _nb_locus, and BinaryStorageBuffer::store().

◆ write_effects_to_parameter()

void TProtoDeletMutations_bitstring::write_effects_to_parameter ( )
310{
311 ostringstream r0, r1;
312
313 for(unsigned int i = 0; i < _nb_locus -1; ++i){
314 // record the effects
315 r0<<_effects[0][i]<<",";
316 r1<<_effects[1][i]<<",";
317 }
318
319 r0<<_effects[0][_nb_locus -1];
320 r1<<_effects[1][_nb_locus -1];
321
322 // save the effects in a parameter; will be added to the .log simulation file
323 string arg = "{{" + r0.str() + "} {" + r1.str() + "}}";
324 get_parameter("delet_effects")->setArg(arg);
325 get_parameter("delet_effects")->setIsSet(true);
326
327}
void setArg(string value)
Sets the parameter's argument.
Definition: param.h:118
void setIsSet(bool value)
Sets the _isSet flag.
Definition: param.h:122

References _effects, _nb_locus, SimComponent::get_parameter(), Param::setArg(), and Param::setIsSet().

Referenced by resetParameterFromSource(), and set_effects().

+ Here is the caller graph for this function:

Member Data Documentation

◆ _back_mutation_rate

double TProtoDeletMutations_bitstring::_back_mutation_rate
private

◆ _continuous_effects

bool TProtoDeletMutations_bitstring::_continuous_effects
private

◆ _dist_p1

double TProtoDeletMutations_bitstring::_dist_p1
private

◆ _dist_p2

double TProtoDeletMutations_bitstring::_dist_p2
private

◆ _dominance

double TProtoDeletMutations_bitstring::_dominance
private

◆ _dominance_model

int TProtoDeletMutations_bitstring::_dominance_model
private

◆ _effects

◆ _fitness_model

unsigned int TProtoDeletMutations_bitstring::_fitness_model
private

Referenced by hatch(), and setParameters().

◆ _fitness_scaling_factor

double TProtoDeletMutations_bitstring::_fitness_scaling_factor
private

◆ _inherit_func_ptr

void(TProtoDeletMutations_bitstring::* TProtoDeletMutations_bitstring::_inherit_func_ptr) (sex_t, bitstring *, bitstring **)
private

Referenced by hatch(), and setParameters().

◆ _init_freq

double TProtoDeletMutations_bitstring::_init_freq
private

Referenced by get_init_freq(), and setParameters().

◆ _mut_rate

double TProtoDeletMutations_bitstring::_mut_rate
private

Referenced by get_mut_rate(), hatch(), and setParameters().

◆ _mutation_model

unsigned int TProtoDeletMutations_bitstring::_mutation_model
private

Referenced by hatch(), and setParameters().

◆ _nb_locus

◆ _reader

TTDeletMutBitstrFH* TProtoDeletMutations_bitstring::_reader
private

◆ _set_effects_func

double(TProtoDeletMutations_bitstring::* TProtoDeletMutations_bitstring::_set_effects_func) (void)
private

◆ _stats

TTDeletMutBitstrSH* TProtoDeletMutations_bitstring::_stats
private

◆ _strength

double TProtoDeletMutations_bitstring::_strength
private

◆ _viability_func_ptr

double(TTDeletMutations_bitstring::* TProtoDeletMutations_bitstring::_viability_func_ptr) (void)
private

◆ _writer

TTDeletMutBitstrFH* TProtoDeletMutations_bitstring::_writer
private

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