[Avida-SVN] r3346 - in development/source: actions cpu main
dk at myxo.css.msu.edu
dk at myxo.css.msu.edu
Mon Jun 29 11:02:10 PDT 2009
Author: dk
Date: 2009-06-29 14:02:06 -0400 (Mon, 29 Jun 2009)
New Revision: 3346
Modified:
development/source/actions/PopulationActions.cc
development/source/actions/PrintActions.cc
development/source/cpu/cHardwareBase.cc
development/source/cpu/cHardwareBase.h
development/source/cpu/cHardwareCPU.cc
development/source/main/cAvidaConfig.h
development/source/main/cEnvironment.cc
development/source/main/cEnvironment.h
development/source/main/cGenome.cc
development/source/main/cGenome.h
development/source/main/cOrganism.cc
development/source/main/cOrganism.h
development/source/main/cPopulation.cc
development/source/main/cPopulation.h
development/source/main/cPopulationCell.cc
development/source/main/cPopulationCell.h
development/source/main/cResource.cc
development/source/main/cResource.h
development/source/main/cStats.cc
development/source/main/cStats.h
Log:
Added rough draft of horizontal gene transfer framework - includes resource, genome metabolism, and insertion.
Modified: development/source/actions/PopulationActions.cc
===================================================================
--- development/source/actions/PopulationActions.cc 2009-06-26 15:12:33 UTC (rev 3345)
+++ development/source/actions/PopulationActions.cc 2009-06-29 18:02:06 UTC (rev 3346)
@@ -3297,6 +3297,33 @@
};
+/*! Diffuse HGT genome fragments.
+
+ If HGT is enabled, each cell will gradually (ok, ok - it'll happen quickly) build up
+ a buffer of genome fragments. This event triggers the diffusion of those fragments.
+
+ For right now, we evaluate each cell in order; fix this if we start seeing diffusion
+ artifacts.
+
+ The actual code for performing the diffusion lives over in cPopulationCell::DiffuseGenomeFragments().
+ */
+class cActionDiffuseHGTGenomeFragments : public cAction {
+public:
+ static const cString GetDescription() { return "Arguments: <none>"; }
+
+ //! Constructor.
+ cActionDiffuseHGTGenomeFragments(cWorld* world, const cString& args) : cAction(world, args) {
+ }
+
+ //! Process this event.
+ void Process(cAvidaContext& ctx) {
+ for(int i=0; i<m_world->GetPopulation().GetSize(); ++i) {
+ m_world->GetPopulation().GetCell(i).DiffuseGenomeFragments();
+ }
+ }
+};
+
+
void RegisterPopulationActions(cActionLibrary* action_lib)
{
action_lib->Register<cActionInject>("Inject");
@@ -3380,6 +3407,7 @@
action_lib->Register<cActionKillNBelowResourceThreshold>("KillNBelowResourceThreshold");
action_lib->Register<cActionKillNAboveResourceThreshold>("KillNAboveResourceThreshold");
+ action_lib->Register<cActionDiffuseHGTGenomeFragments>("DiffuseHGTGenomeFragments");
// @DMB - The following actions are DEPRECATED aliases - These will be removed in 2.7.
action_lib->Register<cActionInject>("inject");
Modified: development/source/actions/PrintActions.cc
===================================================================
--- development/source/actions/PrintActions.cc 2009-06-26 15:12:33 UTC (rev 3345)
+++ development/source/actions/PrintActions.cc 2009-06-29 18:02:06 UTC (rev 3346)
@@ -151,7 +151,10 @@
STATS_OUT_FILE(PrintGroupsFormedData, groupformation.dat);
STATS_OUT_FILE(PrintGroupIds, groupids.dat);
+// hgt information
+STATS_OUT_FILE(PrintHGTData, hgt.dat);
+
#define POP_OUT_FILE(METHOD, DEFAULT) /* 1 */ \
class cAction ## METHOD : public cAction { /* 2 */ \
private: /* 3 */ \
@@ -3176,6 +3179,9 @@
// Group Formation
action_lib->Register<cActionPrintGroupsFormedData>("PrintGroupsFormedData");
action_lib->Register<cActionPrintGroupIds>("PrintGroupIds");
+
+ // hgt
+ action_lib->Register<cActionPrintHGTData>("PrintHGTData");
action_lib->Register<cActionSetVerbose>("VERBOSE");
}
Modified: development/source/cpu/cHardwareBase.cc
===================================================================
--- development/source/cpu/cHardwareBase.cc 2009-06-26 15:12:33 UTC (rev 3345)
+++ development/source/cpu/cHardwareBase.cc 2009-06-29 18:02:06 UTC (rev 3346)
@@ -1077,3 +1077,22 @@
{
m_world->GetDriver().RaiseFatalException(1, "Method cHardwareBase::ReceiveFlash must be overriden.");
}
+
+/*! Retrieve a fragment of this organism's genome that extends downstream from the read head.
+ */
+cGenome cHardwareBase::GetGenomeFragment(unsigned int downstream) {
+ cHeadCPU tmp(GetHead(nHardware::HEAD_READ));
+ cGenome fragment(downstream);
+ for(; downstream>0; --downstream, tmp.Advance()) {
+ fragment.Append(tmp.GetInst());
+ }
+ return fragment;
+}
+
+/*! Insert a genome fragment at the current write head.
+ */
+void cHardwareBase::InsertGenomeFragment(const cGenome& fragment) {
+ cHeadCPU& wh = GetHead(nHardware::HEAD_WRITE);
+ wh.GetMemory().Insert(wh.GetPosition(), fragment);
+ wh.Adjust();
+}
Modified: development/source/cpu/cHardwareBase.h
===================================================================
--- development/source/cpu/cHardwareBase.h 2009-06-26 15:12:33 UTC (rev 3345)
+++ development/source/cpu/cHardwareBase.h 2009-06-29 18:02:06 UTC (rev 3346)
@@ -39,13 +39,13 @@
#ifndef tSmartArray_h
#include "tSmartArray.h"
#endif
+#include "cGenome.h"
using namespace std;
class cAvidaContext;
class cCodeLabel;
class cCPUMemory;
-class cGenome;
class cHardwareTracer;
class cHeadCPU;
class cInjectGenotype;
@@ -193,7 +193,12 @@
// -------- Synchronization --------
//! Called when the organism that owns this CPU has received a flash from a neighbor.
virtual void ReceiveFlash();
-
+
+ // -------- HGT --------
+ //! Retrieve a genome fragment extending downstream from the read head.
+ virtual cGenome GetGenomeFragment(unsigned int downstream);
+ //! Insert a genome fragment at the current write head.
+ virtual void InsertGenomeFragment(const cGenome& fragment);
protected:
// -------- Core Execution Methods --------
Modified: development/source/cpu/cHardwareCPU.cc
===================================================================
--- development/source/cpu/cHardwareCPU.cc 2009-06-26 15:12:33 UTC (rev 3345)
+++ development/source/cpu/cHardwareCPU.cc 2009-06-29 18:02:06 UTC (rev 3346)
@@ -5746,7 +5746,8 @@
} else
doSlipMutation(ctx, write_head.GetMemory(), write_head.GetPosition());
}
-
+ m_organism->AttemptHGTInsertion(ctx);
+
read_head.Advance();
write_head.Advance();
return true;
Modified: development/source/main/cAvidaConfig.h
===================================================================
--- development/source/main/cAvidaConfig.h 2009-06-26 15:12:33 UTC (rev 3345)
+++ development/source/main/cAvidaConfig.h 2009-06-29 18:02:06 UTC (rev 3346)
@@ -647,10 +647,21 @@
CONFIG_ADD_VAR(USE_FORM_GROUPS, int, 0, "Enable organisms to form groups. 0=off,\n 1=on no restrict,\n 2=on restrict to defined");
// -------- Deme network config options --------
+ CONFIG_ADD_GROUP(DEME_NETWORK_GROUP, "Deme network settings");
CONFIG_ADD_VAR(DEME_NETWORK_TYPE, int, 0, "0=topology, structure of network determines fitness.");
CONFIG_ADD_VAR(DEME_NETWORK_REQUIRES_CONNECTEDNESS, int, 1, "Whether the deme's network must be connected before an actual fitness is calculated.");
CONFIG_ADD_VAR(DEME_NETWORK_TOPOLOGY_FITNESS, int, 0, "Network measure used to determine fitness; see cDemeTopologyNetwork.h.");
+ // -------- Horizontal gene transfer config options --------
+ CONFIG_ADD_GROUP(HGT_GROUP, "Horizontal gene transfer settings");
+ CONFIG_ADD_VAR(ENABLE_HGT, int, 0, "Whether HGT is enabled; 0=false (default),\n 1=true.");
+ CONFIG_ADD_VAR(HGT_FRAGMENT_SIZE_MEAN, double, 10, "Mean size of fragments (drawn from a normal\ndist., default=10).");
+ CONFIG_ADD_VAR(HGT_FRAGMENT_SIZE_VARIANCE, double, 2, "Variance of fragments (drawn from a normal\ndist., default=2).");
+ CONFIG_ADD_VAR(HGT_MAX_FRAGMENTS_PER_CELL, int, 100, "Max. allowed number of fragments\nper cell (default=100).");
+ CONFIG_ADD_VAR(HGT_DIFFUSION_METHOD, int, 0, "Method to use for diffusion of genome\nfragments (0=none [default]).");
+ CONFIG_ADD_VAR(HGT_INSERTION_PROB, double, 0.0, "Probability that a genome fragment\nwill be inserted during a copy (default=0.0).");
+ CONFIG_ADD_VAR(HGT_LOOKAHEAD_LENGTH, int, 4, "Number of instructions forward from the\nread head that will be used to calculate\nthe liklihood of HGT.");
+
#endif
inline void Load(const cString& filename) { Load(filename, false); }
Modified: development/source/main/cEnvironment.cc
===================================================================
--- development/source/main/cEnvironment.cc 2009-06-26 15:12:33 UTC (rev 3345)
+++ development/source/main/cEnvironment.cc 2009-06-29 18:02:06 UTC (rev 3346)
@@ -29,6 +29,7 @@
#include "cEnvironment.h"
+#include "defs.h"
#include "cArgSchema.h"
#include "cAvidaContext.h"
#include "cEnvReqs.h"
@@ -36,6 +37,8 @@
#include "cInitFile.h"
#include "cInstSet.h"
#include "nMutation.h"
+#include "cPopulation.h"
+#include "cPopulationCell.h"
#include "cRandom.h"
#include "cReaction.h"
#include "nReaction.h"
@@ -63,6 +66,13 @@
using namespace std;
+cEnvironment::cEnvironment(cWorld* world) : m_world(world) , m_tasklib(world),
+m_input_size(INPUT_SIZE_DEFAULT), m_output_size(OUTPUT_SIZE_DEFAULT), m_true_rand(false),
+m_use_specific_inputs(false), m_specific_inputs(), m_mask(0)
+{
+ mut_rates.Setup(world);
+}
+
cEnvironment::~cEnvironment()
{
for (int i = 0; i < m_state_grids.GetSize(); i++) delete m_state_grids[i];
@@ -452,13 +462,41 @@
cerr <<"Error: Energy resources can not be used without the energy model.\n";
}
}
+ else if (var_name == "hgt") {
+ // this resource is for HGT -- corresponds to genome fragments present in cells.
+ if(!AssertInputBool(var_value, "hgt", var_type)) return false;
+ new_resource->SetHGTMetabolize(var_value.AsInt());
+ }
else {
cerr << "Error: Unknown variable '" << var_name
<< "' in resource '" << name << "'" << endl;
return false;
}
}
-
+
+ // Prevent misconfiguration of HGT:
+ if(new_resource->GetHGTMetabolize() &&
+ ((new_resource->GetGeometry() != nGeometry::GLOBAL)
+ || (new_resource->GetInitial() > 0.0)
+ || (new_resource->GetInflow() > 0.0)
+ || (new_resource->GetOutflow() > 0.0)
+ || (new_resource->GetInflowX1() != -99)
+ || (new_resource->GetInflowX2() != -99)
+ || (new_resource->GetInflowY1() != -99)
+ || (new_resource->GetInflowY2() != -99)
+ || (new_resource->GetXDiffuse() != 1.0)
+ || (new_resource->GetXGravity() != 0.0)
+ || (new_resource->GetYDiffuse() != 1.0)
+ || (new_resource->GetYGravity() != 0.0)
+ || (new_resource->GetDemeResource() != false))) {
+ cerr << "Error: misconfigured HGT resource: " << name << endl;
+ return false;
+ }
+ if(new_resource->GetHGTMetabolize() && !m_world->GetConfig().ENABLE_HGT.Get()) {
+ cerr << "Error: resource configured to use HGT, but HGT not enabled." << endl;
+ return false;
+ }
+
// If there are valid values for X/Y1's but not for X/Y2's assume that
// the user is interested only in one point and set the X/Y2's to the
// same value as X/Y1's
@@ -1077,7 +1115,7 @@
// And let's process it!
DoProcesses(ctx, cur_reaction->GetProcesses(), resource_count, rbins_count,
- task_quality, task_probability, task_cnt, i, result);
+ task_quality, task_probability, task_cnt, i, result, taskctx);
// Note: the reaction is actually marked as being performed inside DoProcesses.
@@ -1182,7 +1220,7 @@
void cEnvironment::DoProcesses(cAvidaContext& ctx, const tList<cReactionProcess>& process_list,
const tArray<double>& resource_count, const tArray<double>& rbins_count,
const double task_quality, const double task_probability, const int task_count,
- const int reaction_id, cReactionResult& result) const
+ const int reaction_id, cReactionResult& result, cTaskContext& taskctx) const
{
const int num_process = process_list.GetSize();
@@ -1209,7 +1247,31 @@
if (in_resource == NULL) {
// Test if infinite resource
consumed = max_consumed * local_task_quality * task_plasticity_modifier;
- } else {
+
+ } else if(in_resource->GetHGTMetabolize()) {
+ /* HGT Metabolism
+ This bit of code is triggered when ENABLE_HGT=1 and a resource has hgt=1.
+ Here's the idea: Each cell in the environment holds a buffer of genome fragments,
+ where these fragments are drawn from the remains of organisms that have died.
+ These remains are a potential source of energy to the current inhabitant of the
+ cell. This code metabolizes one of those fragments by pretending that it's just
+ another resource. Task quality can be used to control the conversion of fragments
+ to bonus, but the amount of resource consumed is always equal to the length of the
+ fragment.
+ */
+ int cellid = taskctx.GetOrganism()->GetCellID();
+ if(cellid != -1) { // can't do this in the test cpu
+ cPopulationCell& cell = m_world->GetPopulation().GetCell(cellid);
+ if(cell.CountGenomeFragments() > 0) {
+ cGenome fragment = cell.PopGenomeFragment();
+ consumed = local_task_quality * fragment.GetSize();
+ result.Consume(in_resource->GetID(), fragment.GetSize(), true);
+ m_world->GetStats().GenomeFragmentMetabolized(taskctx.GetOrganism(), fragment);
+ }
+ }
+ // if we can't metabolize a fragment, stop here.
+ if(consumed == 0.0) { continue; }
+ } else {
// Otherwise we're using a finite resource
const int res_id = in_resource->GetID();
@@ -1265,7 +1327,7 @@
result.Consume(res_id, consumed, !using_rbins);
}
}
-
+
// Mark the reaction as having been performed if we get here.
result.MarkReaction(reaction_id);
Modified: development/source/main/cEnvironment.h
===================================================================
--- development/source/main/cEnvironment.h 2009-06-26 15:12:33 UTC (rev 3345)
+++ development/source/main/cEnvironment.h 2009-06-29 18:02:06 UTC (rev 3346)
@@ -48,9 +48,6 @@
#ifndef cTaskLib_h
#include "cTaskLib.h"
#endif
-#ifndef defs_h
-#include "defs.h"
-#endif
#ifndef tList_h
#include "tList.h"
#endif
@@ -127,14 +124,14 @@
const tArray<double>& resource_count, const tArray<double>& rbin_count,
const double task_quality, const double task_probability,
const int task_count, const int reaction_id,
- cReactionResult& result) const;
+ cReactionResult& result, cTaskContext& taskctx) const;
cEnvironment(); // @not_implemented
cEnvironment(const cEnvironment&); // @not_implemented
cEnvironment& operator=(const cEnvironment&); // @not_implemented
public:
- inline cEnvironment(cWorld* world);
+ cEnvironment(cWorld* world);
~cEnvironment();
bool Load(const cString& filename); // Reads the environment from disk.
@@ -202,15 +199,6 @@
};
-inline cEnvironment::cEnvironment(cWorld* world) : m_world(world) , m_tasklib(world),
- m_input_size(INPUT_SIZE_DEFAULT), m_output_size(OUTPUT_SIZE_DEFAULT), m_true_rand(false),
- m_use_specific_inputs(false), m_specific_inputs(), m_mask(0)
-{
- mut_rates.Setup(world);
-}
-
-
-
#ifdef ENABLE_UNIT_TESTS
namespace nEnvironment {
/**
Modified: development/source/main/cGenome.cc
===================================================================
--- development/source/main/cGenome.cc 2009-06-26 15:12:33 UTC (rev 3345)
+++ development/source/main/cGenome.cc 2009-06-29 18:02:06 UTC (rev 3346)
@@ -24,6 +24,7 @@
*/
#include "cGenome.h"
+#include <iterator>
using namespace std;
@@ -62,11 +63,12 @@
It expects STL semantics for an iterator range. We're avoiding templating this
(for now). Refactor if a new range type is needed.
*/
-cGenome::cGenome(cInstruction* begin, cInstruction* end) : m_active_size(0)
+cGenome::cGenome(const cInstruction* begin, const cInstruction* end) : m_active_size(0)
{
- m_genome.Resize((end - begin) / sizeof(cInstruction));
-
- for (cInstruction* i = begin; i != end; i++, m_active_size++) m_genome[m_active_size] = *i;
+ m_genome.Resize(std::distance(begin,end));
+ for(const cInstruction* i=begin; i!=end; i++, m_active_size++) {
+ m_genome[m_active_size] = *i;
+ }
}
Modified: development/source/main/cGenome.h
===================================================================
--- development/source/main/cGenome.h 2009-06-26 15:12:33 UTC (rev 3345)
+++ development/source/main/cGenome.h 2009-06-29 18:02:06 UTC (rev 3346)
@@ -63,7 +63,7 @@
explicit cGenome(int _size); //! Constructor that builds a 'blank' cGenome of the specified size
cGenome(const cGenome& in_genome); //! Copy constructor
cGenome(const cString& in_string); //! Constructor that builds genome from a string
- cGenome(cInstruction* begin, cInstruction* end); //! Constructor that builds genome from a range of instructions
+ cGenome(const cInstruction* begin, const cInstruction* end); //! Constructor that builds genome from a range of instructions
virtual ~cGenome(); //! Virtual destructor; there are subclasses.
inline int GetSize() const { return m_active_size; }
@@ -89,7 +89,6 @@
virtual bool operator==(const cGenome& other_genome) const;
virtual bool operator!=(const cGenome& other_genome) const { return !(this->operator==(other_genome)); }
virtual bool operator<(const cGenome& other_genome) const { return AsString() < other_genome.AsString(); }
-
};
#endif
Modified: development/source/main/cOrganism.cc
===================================================================
--- development/source/main/cOrganism.cc 2009-06-26 15:12:33 UTC (rev 3345)
+++ development/source/main/cOrganism.cc 2009-06-29 18:02:06 UTC (rev 3346)
@@ -1234,3 +1234,48 @@
return val;
}
+
+/*! Tests for and attempts to perform an insertion of an HGT genome fragment
+ into this organism's genome. Returns true if a genome fragment was inserted,
+ false otherwise.
+ */
+bool cOrganism::AttemptHGTInsertion(cAvidaContext& ctx) {
+ if(!m_world->GetConfig().ENABLE_HGT.Get()
+ || (m_world->GetConfig().HGT_INSERTION_PROB.Get()==0.0)) {
+ return false;
+ }
+
+ int cellid = GetCellID();
+ if(cellid == -1) { return false; } // test cpu, nothing to do.
+ cPopulationCell& cell = m_world->GetPopulation().GetCell(cellid);
+
+ // need to calculate the probability of inserting any of the fragments that live in
+ // this organism's cell. this P() is partially based on the instructions near the
+ // read head, so let's grab that fragment:
+ cGenome nearby = m_hardware->GetGenomeFragment(m_world->GetConfig().HGT_LOOKAHEAD_LENGTH.Get());
+
+ // get the list of genome fragments currently in this cell:
+ typedef cPopulationCell::fragment_list_type fragment_list_type;
+ fragment_list_type& fragments = cell.GetFragments();
+
+ // now, the probability of selecting a fragment for insertion is a combination
+ // of how long it's been in contact with this organism, a background rate, and
+ // the similarity of this fragment to the nearby fragment. the oldest genomes
+ // have the greatest chance of being inserted, because they've gone through this
+ // loop multiple times.
+ //
+ // when a fragment is inserted, we also need to adjust the amount of hgt resource
+ // present in the population, and do some stat tracking.
+ for(fragment_list_type::iterator i=fragments.begin(); i!=fragments.end(); ++i) {
+ int d = std::max(cGenomeUtil::FindHammingDistance(nearby, *i), 1);
+ if(m_world->GetRandom().P(m_world->GetConfig().HGT_INSERTION_PROB.Get() / d)) {
+ m_hardware->InsertGenomeFragment(*i);
+ m_world->GetPopulation().AdjustHGTResource(-i->GetSize());
+ m_world->GetStats().GenomeFragmentInserted(this, *i);
+ fragments.erase(i);
+ return true;
+ }
+ }
+
+ return false;
+}
Modified: development/source/main/cOrganism.h
===================================================================
--- development/source/main/cOrganism.h 2009-06-26 15:12:33 UTC (rev 3345)
+++ development/source/main/cOrganism.h 2009-06-29 18:02:06 UTC (rev 3346)
@@ -321,7 +321,6 @@
cInjectGenotype& GetParasite(int x) { return *m_parasites[x]; }
int GetNumParasites() const { return m_parasites.GetSize(); }
void ClearParasites();
-
// -------- Mutation Rate Convenience Methods --------
bool TestCopyMut(cAvidaContext& ctx) const { return m_mut_rates.TestCopyMut(ctx); }
@@ -329,7 +328,7 @@
bool TestCopyDel(cAvidaContext& ctx) const { return m_mut_rates.TestCopyDel(ctx); }
bool TestCopyUniform(cAvidaContext& ctx) const { return m_mut_rates.TestCopyUniform(ctx); }
bool TestCopySlip(cAvidaContext& ctx) const { return m_mut_rates.TestCopySlip(ctx); }
-
+
bool TestDivideMut(cAvidaContext& ctx) const { return m_mut_rates.TestDivideMut(ctx); }
bool TestDivideIns(cAvidaContext& ctx) const { return m_mut_rates.TestDivideIns(ctx); }
bool TestDivideDel(cAvidaContext& ctx) const { return m_mut_rates.TestDivideDel(ctx); }
@@ -647,6 +646,11 @@
/*! The main DoOutput function. The DoOutputs above all forward to this function. */
void doOutput(cAvidaContext& ctx, tBuffer<int>& input_buffer, tBuffer<int>& output_buffer, const bool on_divide);
+
+ // -------- HGT methods --------
+public:
+ //! Tests for and attempts to perform an insertion of an HGT genome fragment into this organism's genome.
+ bool AttemptHGTInsertion(cAvidaContext& ctx);
};
Modified: development/source/main/cPopulation.cc
===================================================================
--- development/source/main/cPopulation.cc 2009-06-26 15:12:33 UTC (rev 3345)
+++ development/source/main/cPopulation.cc 2009-06-29 18:02:06 UTC (rev 3346)
@@ -91,6 +91,7 @@
, environment(world->GetEnvironment())
, num_organisms(0)
, sync_events(false)
+, m_hgt_resid(-1)
{
// Avida specific information.
world_x = world->GetConfig().WORLD_X.Get();
@@ -243,6 +244,15 @@
for (int i = 0; i < resource_lib.GetSize(); i++) {
cResource * res = resource_lib.GetResource(i);
+
+ // check to see if this is the hgt resource:
+ if(res->GetHGTMetabolize()) {
+ if(m_hgt_resid != -1) {
+ m_world->GetDriver().RaiseFatalException(-1, "Only one HGT resource is currently supported.");
+ }
+ m_hgt_resid = i;
+ }
+
if (!res->GetDemeResource()) {
global_res_index++;
const double decay = 1.0 - res->GetOutflow();
@@ -268,6 +278,11 @@
exit(1);
}
}
+
+ // if HGT is on, make sure there's a resource for it:
+ if(m_world->GetConfig().ENABLE_HGT.Get() && (m_hgt_resid == -1)) {
+ m_world->GetDriver().RaiseFatalException(-1, "HGT is enabled, but no HGT resource is defined; add hgt=1 to a single resource in the environment file.");
+ }
}
@@ -803,8 +818,14 @@
deme_array[in_cell.GetDemeID()].OrganismDeath(in_cell);
}
genotype->RemoveOrganism();
-
organism->ClearParasites();
+
+ // If HGT is turned on, this organism's genome needs to be split up into fragments
+ // and deposited in its cell. We then also have to add the size of this genome to
+ // the HGT resource.
+ if(m_world->GetConfig().ENABLE_HGT.Get()) {
+ in_cell.AddGenomeFragments(organism->GetGenome());
+ }
// And clear it!
in_cell.RemoveOrganism();
@@ -5699,3 +5720,9 @@
return num_orgs;
}
+
+/*! Modify current level of the HGT resource.
+ */
+void cPopulation::AdjustHGTResource(double delta) {
+ resource_count.Modify(m_hgt_resid, delta);
+}
Modified: development/source/main/cPopulation.h
===================================================================
--- development/source/main/cPopulation.h 2009-06-26 15:12:33 UTC (rev 3345)
+++ development/source/main/cPopulation.h 2009-06-29 18:02:06 UTC (rev 3346)
@@ -118,7 +118,6 @@
// Group formation information
std::map<int, int> m_groups; //<! Maps the group id to the number of orgs in the group
-
///////////////// Private Methods ////////////////////
void BuildTimeSlicer(cChangeList* change_list); // Build the schedule object
@@ -366,6 +365,12 @@
inline bool operator>=(const sTmpGenotype& rhs) const { return id_num >= rhs.id_num; }
};
+ // -------- HGT support --------
+private:
+ int m_hgt_resid; //!< HGT resource ID.
+public:
+ //! Modify current level of the HGT resource.
+ void AdjustHGTResource(double delta);
};
Modified: development/source/main/cPopulationCell.cc
===================================================================
--- development/source/main/cPopulationCell.cc 2009-06-26 15:12:33 UTC (rev 3345)
+++ development/source/main/cPopulationCell.cc 2009-06-29 18:02:06 UTC (rev 3346)
@@ -25,6 +25,7 @@
#include "cPopulationCell.h"
+#include "cDoubleSum.h"
#include "nHardware.h"
#include "cOrganism.h"
#include "cTools.h"
@@ -34,6 +35,9 @@
#include "cPopulation.h"
#include "cDeme.h"
+#include <cmath>
+#include <iterator>
+
using namespace std;
@@ -46,6 +50,7 @@
, m_deme_id(in_cell.m_deme_id)
, m_cell_data(in_cell.m_cell_data)
, m_spec_state(in_cell.m_spec_state)
+ , m_hgt(0)
{
// Copy the mutation rates into a new structure
m_mut_rates = new cMutationRates(*in_cell.m_mut_rates);
@@ -54,29 +59,45 @@
tConstListIterator<cPopulationCell> conn_it(in_cell.m_connections);
cPopulationCell* test_cell;
while ((test_cell = const_cast<cPopulationCell*>(conn_it.Next()))) m_connections.PushRear(test_cell);
+
+ // copy the hgt information, if needed.
+ if(in_cell.m_hgt) {
+ InitHGTSupport();
+ *m_hgt = *in_cell.m_hgt;
+ }
}
void cPopulationCell::operator=(const cPopulationCell& in_cell)
{
- m_world = in_cell.m_world;
- m_organism = in_cell.m_organism;
- m_hardware = in_cell.m_hardware;
- m_inputs = in_cell.m_inputs;
- m_cell_id = in_cell.m_cell_id;
- m_deme_id = in_cell.m_deme_id;
- m_cell_data = in_cell.m_cell_data;
- m_spec_state = in_cell.m_spec_state;
+ if(this != &in_cell) {
+ m_world = in_cell.m_world;
+ m_organism = in_cell.m_organism;
+ m_hardware = in_cell.m_hardware;
+ m_inputs = in_cell.m_inputs;
+ m_cell_id = in_cell.m_cell_id;
+ m_deme_id = in_cell.m_deme_id;
+ m_cell_data = in_cell.m_cell_data;
+ m_spec_state = in_cell.m_spec_state;
- // Copy the mutation rates, constructing the structure as necessary
- if (m_mut_rates == NULL)
- m_mut_rates = new cMutationRates(*in_cell.m_mut_rates);
- else
- m_mut_rates->Copy(*in_cell.m_mut_rates);
+ // Copy the mutation rates, constructing the structure as necessary
+ if (m_mut_rates == NULL)
+ m_mut_rates = new cMutationRates(*in_cell.m_mut_rates);
+ else
+ m_mut_rates->Copy(*in_cell.m_mut_rates);
- // Copy the connection list
- tConstListIterator<cPopulationCell> conn_it(in_cell.m_connections);
- cPopulationCell* test_cell;
- while ((test_cell = const_cast<cPopulationCell*>(conn_it.Next()))) m_connections.PushRear(test_cell);
+ // Copy the connection list
+ tConstListIterator<cPopulationCell> conn_it(in_cell.m_connections);
+ cPopulationCell* test_cell;
+ while ((test_cell = const_cast<cPopulationCell*>(conn_it.Next()))) m_connections.PushRear(test_cell);
+
+ // copy hgt information, if needed.
+ delete m_hgt;
+ m_hgt = 0;
+ if(in_cell.m_hgt) {
+ InitHGTSupport();
+ *m_hgt = *in_cell.m_hgt;
+ }
+ }
}
void cPopulationCell::Setup(cWorld* world, int in_id, const cMutationRates& in_rates, int x, int y)
@@ -255,3 +276,87 @@
// Nothing for the moment...
return true;
}
+
+/*! Diffuse genome fragments from this cell to its neighbors.
+
+ NOTE: This method is for OUTGOING diffusion only.
+
+ There are many possible ways in which genome fragments could be diffused. We'll
+ put in the framework to support those other mechanisms, but we're not going to
+ worry about this until we need it. Not terribly interested in recreating an
+ artificial chemistry here...
+ */
+void cPopulationCell::DiffuseGenomeFragments() {
+ InitHGTSupport();
+
+ switch(m_world->GetConfig().HGT_DIFFUSION_METHOD.Get()) {
+ case 0: { // none
+ break;
+ }
+ default: {
+ m_world->GetDriver().RaiseFatalException(-1, "Unrecognized diffusion type in cPopulationCell::DiffuseGenomeFragments().");
+ }
+ }
+}
+
+/*! Add fragments from the passed-in genome to the HGT fragments contained in this cell.
+
+ Split the passed-in genome into fragments according to a normal distribution specified
+ by HGT_FRAGMENT_SIZE_MEAN and HGT_FRAGMENT_SIZE_VARIANCE. These fragments are added
+ to this cell's fragment list.
+
+ As a safety measure, we also remove old fragments to conserve memory. Specifically, we
+ remove old fragments until at most HGT_MAX_FRAGMENTS_PER_CELL fragments remain.
+ */
+void cPopulationCell::AddGenomeFragments(const cGenome& genome) {
+ assert(genome.GetSize()>0); // oh, sweet sanity.
+ InitHGTSupport();
+
+ m_world->GetPopulation().AdjustHGTResource(genome.GetSize());
+
+ // chop this genome up into pieces, add each to the back of this cell's buffer.
+ int remaining_size=genome.GetSize();
+ const cInstruction* i=&genome[0];
+ do {
+ int fsize = std::min(remaining_size,
+ (int)floor(m_world->GetRandom().GetRandNormal(m_world->GetConfig().HGT_FRAGMENT_SIZE_MEAN.Get(),
+ m_world->GetConfig().HGT_FRAGMENT_SIZE_VARIANCE.Get())));
+ m_hgt->fragments.push_back(cGenome(i, i+fsize));
+ i+=fsize;
+ remaining_size-=fsize;
+ } while(remaining_size>0);
+
+ // pop off the front of this cell's buffer until we have <= HGT_MAX_FRAGMENTS_PER_CELL.
+ while(m_hgt->fragments.size()>(unsigned int)m_world->GetConfig().HGT_MAX_FRAGMENTS_PER_CELL.Get()) {
+ m_world->GetPopulation().AdjustHGTResource(-m_hgt->fragments.front().GetSize());
+ m_hgt->fragments.pop_front();
+ }
+}
+
+/*! Retrieve the number of genome fragments currently found in this cell.
+ */
+unsigned int cPopulationCell::CountGenomeFragments() const {
+ if(IsHGTInitialized()) {
+ return m_hgt->fragments.size();
+ } else {
+ return 0;
+ }
+}
+
+/*! Remove and return a random genome fragment.
+ */
+cGenome cPopulationCell::PopGenomeFragment() {
+ assert(m_hgt!=0);
+ fragment_list_type::iterator i = m_hgt->fragments.begin();
+ std::advance(i, m_world->GetRandom().GetUInt(0, m_hgt->fragments.size()));
+ cGenome tmp = *i;
+ m_hgt->fragments.erase(i);
+ return tmp;
+}
+
+/*! Retrieve the list of fragments from this cell.
+ */
+cPopulationCell::fragment_list_type& cPopulationCell::GetFragments() {
+ InitHGTSupport();
+ return m_hgt->fragments;
+}
Modified: development/source/main/cPopulationCell.h
===================================================================
--- development/source/main/cPopulationCell.h 2009-06-26 15:12:33 UTC (rev 3345)
+++ development/source/main/cPopulationCell.h 2009-06-29 18:02:06 UTC (rev 3346)
@@ -28,6 +28,7 @@
#include <fstream>
#include <set>
+#include <deque>
#ifndef cMutationRates_h
#include "cMutationRates.h"
@@ -38,6 +39,7 @@
#ifndef tList_h
#include "tList.h"
#endif
+#include "cGenome.h"
class cHardwareBase;
class cPopulation;
@@ -78,9 +80,11 @@
public:
- cPopulationCell() : m_world(NULL), m_organism(NULL), m_hardware(NULL), m_mut_rates(NULL), m_migrant(false) { ; }
+ typedef std::set<cPopulationCell*> neighborhood_type; //!< Type for cell neighborhoods.
+
+ cPopulationCell() : m_world(NULL), m_organism(NULL), m_hardware(NULL), m_mut_rates(NULL), m_migrant(false), m_hgt(0) { ; }
cPopulationCell(const cPopulationCell& in_cell);
- ~cPopulationCell() { delete m_mut_rates; }
+ ~cPopulationCell() { delete m_mut_rates; delete m_hgt; }
void operator=(const cPopulationCell& in_cell);
@@ -126,6 +130,31 @@
double UptakeCellEnergy(double frac_to_uptake);
bool OK();
+
+ // -------- HGT support --------
+public:
+ typedef std::deque<cGenome> fragment_list_type; //!< Type for the list of genome fragments.
+ //! Diffuse genome fragments from this cell to its neighbors.
+ void DiffuseGenomeFragments();
+ //! Add fragments from the passed-in genome to the HGT fragments contained in this cell.
+ void AddGenomeFragments(const cGenome& genome);
+ //! Retrieve the number of genome fragments currently found in this cell.
+ unsigned int CountGenomeFragments() const;
+ //! Remove and return the front genome fragment.
+ cGenome PopGenomeFragment();
+ //! Retrieve the list of fragments from this cell.
+ fragment_list_type& GetFragments();
+private:
+ //! Contains HGT-related information for this cell.
+ struct HGTSupport {
+ // WARNING: the default operator= is used in cPopulationCell's copy ctor and operator=.
+ fragment_list_type fragments; //!< Fragments located in this cell.
+ };
+ HGTSupport* m_hgt; //!< Lazily-initialized pointer to the HGT support struct.
+ //! Initialize HGT support in this cell.
+ inline void InitHGTSupport() { if(!m_hgt) { m_hgt = new HGTSupport(); } }
+ //! Is HGT initialized?
+ inline bool IsHGTInitialized() const { return m_hgt != 0; }
};
inline int cPopulationCell::GetInputAt(int& input_pointer)
@@ -134,6 +163,4 @@
return m_inputs[input_pointer++];
}
-
-
#endif
Modified: development/source/main/cResource.cc
===================================================================
--- development/source/main/cResource.cc 2009-06-26 15:12:33 UTC (rev 3345)
+++ development/source/main/cResource.cc 2009-06-29 18:02:06 UTC (rev 3346)
@@ -67,6 +67,7 @@
, ygravity(0.0)
, deme_resource(false)
, energy_resource(false)
+ , hgt_metabolize(false)
{
}
Modified: development/source/main/cResource.h
===================================================================
--- development/source/main/cResource.h 2009-06-26 15:12:33 UTC (rev 3345)
+++ development/source/main/cResource.h 2009-06-29 18:02:06 UTC (rev 3346)
@@ -85,7 +85,7 @@
bool energy_resource; // only implemented for spacial resource
tArray<cCellResource> cell_list;
tArray<int> cell_id_list;
-
+ bool hgt_metabolize;
cResource(); // @not_implemented
@@ -115,8 +115,8 @@
bool GetEnergyResource() const { return energy_resource; }
tArray<cCellResource> *GetCellListPtr() { return &cell_list; }
tArray<int> *GetCellIdListPtr() { return &cell_id_list; }
+ bool GetHGTMetabolize() const { return hgt_metabolize; }
-
void SetInitial(double _initial) { initial = _initial; }
void SetInflow (double _inflow ) { inflow = _inflow; }
void SetOutflow(double _outflow) { outflow = _outflow; }
@@ -140,6 +140,7 @@
void UpdateCellResource(cCellResource *_CellResoucePtr, double _initial,
double _inflow, double _outflow);
void SetCellIdList(tArray<int>& id_list); //SLG partial resources
+ void SetHGTMetabolize(int _in) { hgt_metabolize = _in; }
};
Modified: development/source/main/cStats.cc
===================================================================
--- development/source/main/cStats.cc 2009-06-26 15:12:33 UTC (rev 3345)
+++ development/source/main/cStats.cc 2009-06-29 18:02:06 UTC (rev 3346)
@@ -2775,3 +2775,31 @@
df.WriteComment("Deme network statistics");
df.WriteTimeStamp();
}
+
+/*! Called when an organism metabolizes a genome fragment.
+ */
+void cStats::GenomeFragmentMetabolized(cOrganism* organism, const cGenome& fragment) {
+ m_hgt_metabolized.Add(fragment.GetSize());
+}
+
+void cStats::GenomeFragmentInserted(cOrganism* organism, const cGenome& fragment) {
+ m_hgt_inserted.Add(fragment.GetSize());
+}
+
+/*! Print HGT statistics.
+ */
+void cStats::PrintHGTData(const cString& filename) {
+ cDataFile& df = m_world->GetDataFile(filename);
+
+ df.WriteComment("Horizontal gene transfer statistics");
+ df.WriteTimeStamp();
+ df.Write(GetUpdate(), "Update [update]");
+ df.Write(m_hgt_metabolized.Count(), "Total count of metabolized genome fragments [metcount]");
+ df.Write(m_hgt_metabolized.Sum(), "Total size of metabolized genome fragments [metsize]");
+ df.Write(m_hgt_inserted.Count(), "Total count of insertion events [inscount]");
+ df.Write(m_hgt_inserted.Sum(), "Total size of insertion events [inssize]");
+ df.Endl();
+
+ m_hgt_metabolized.Clear();
+ m_hgt_inserted.Clear();
+}
Modified: development/source/main/cStats.h
===================================================================
--- development/source/main/cStats.h 2009-06-26 15:12:33 UTC (rev 3345)
+++ development/source/main/cStats.h 2009-06-29 18:02:06 UTC (rev 3346)
@@ -59,6 +59,7 @@
#ifndef nGeometry_h
#include "nGeometry.h"
#endif
+#include "cGenome.h"
#if USE_tMemTrack
# ifndef tMemTrack_h
@@ -972,6 +973,18 @@
}
//! Print network statistics.
void PrintDemeNetworkData(const cString& filename);
+
+ // -------- HGT support --------
+private:
+ cDoubleSum m_hgt_metabolized; //!< Total length of metabolized genome fragments.
+ cDoubleSum m_hgt_inserted; //!< Total length of inserted genome fragments.
+public:
+ //! Called when an organism metabolizes a genome fragment.
+ void GenomeFragmentMetabolized(cOrganism* organism, const cGenome& fragment);
+ //! Called when an organism inserts a genome fragment.
+ void GenomeFragmentInserted(cOrganism* organism, const cGenome& fragment);
+ //! Print HGT statistics.
+ void PrintHGTData(const cString& filename);
};
More information about the Avida-cvs
mailing list