Mentions légales du service

Skip to content
Snippets Groups Projects

Compare revisions

Changes are shown as if the source revision was being merged into the target revision. Learn more about comparing revisions.

Source

Select target project
No results found

Target

Select target project
  • rouzaudc/aevol
  • aevol/aevol
  • tgrohens/aevol
  • mfoley/aevol
  • jluisell/aevol
  • jluisell/aevol-eukaryotes
  • rouzaudc/aevol-asyncgarbage
  • elie.dumont/aevol
8 results
Show changes
Hide whitespace changes
Inline Side-by-side
Showing
with 639 additions and 837 deletions
src/libaevol/tests/Test_ae_list.cpp deleted 100644 → 0
View file @ 5d97e536
// ****************************************************************************
//
// Aevol - An in silico experimental evolution platform
//
// ****************************************************************************
//
// Copyright: See the AUTHORS file provided with the package or <www.aevol.fr>
// Web: http://www.aevol.fr/
// E-mail: See <http://www.aevol.fr/contact/>
// Original Authors : Guillaume Beslon, Carole Knibbe, David Parsons
//
// This program is free software: you can redistribute it and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation, either version 2 of the License, or
// (at your option) any later version.
//
// This program is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
// GNU General Public License for more details.
//
// You should have received a copy of the GNU General Public License
// along with this program. If not, see <http://www.gnu.org/licenses/>.
//
//*****************************************************************************
// =================================================================
// Libraries
// =================================================================
#include <string>
// =================================================================
// Project Files
// =================================================================
#include "Test_ae_list.h"
namespace aevol {
// ===========================================================================
// Declare Used Namespaces
// ===========================================================================
//############################################################################
// #
// Class Test_ae_list #
// #
//############################################################################
CPPUNIT_TEST_SUITE_REGISTRATION(Test_ae_list);
// ===========================================================================
// Static attributes
// ===========================================================================
// Don't change value (hard use in some tests)
const int Test_ae_list::INT_LIST_SIZE = 20;
// ===========================================================================
// Constructors
// ===========================================================================
Test_ae_list::Test_ae_list(void)
{
}
// ===========================================================================
// Destructors
// ===========================================================================
Test_ae_list::~Test_ae_list(void)
{
}
// ===========================================================================
// Operators
// ===========================================================================
// ===========================================================================
// Public Methods
// ===========================================================================
void Test_ae_list::setUp(void)
{
printf("Test_ae_list setUp\n");
int_list = new ae_list<int*>();
for (int i = 0 ; i < INT_LIST_SIZE ; i++)
{
int_list->add(new int(i + 1));
}
}
void Test_ae_list::tearDown(void)
{
int_list->erase(true);
delete int_list;
}
void Test_ae_list::basic_tests1(void)
{
printf("Test_ae_list basic_tests1\n");
// Manually check the content of int_list (checks add(T*))
ae_list_node<int*>* node = int_list->get_first();
int32_t nb_elts = 0;
while (node != NULL)
{
CPPUNIT_ASSERT_EQUAL(++nb_elts, *node->get_obj());
node = node->get_next();
}
CPPUNIT_ASSERT_EQUAL(INT_LIST_SIZE, nb_elts);
// Construct the same list using add_front
ae_list<int*>* expected = new ae_list<int*>();
for (int i = INT_LIST_SIZE ; i > 0 ; i--)
{
expected->add_front(new int(i));
}
assert_equal(expected, int_list, CPPUNIT_SOURCELINE());
expected->erase(true);
delete expected;
// Check erase and is_empty
int_list->erase(true);
CPPUNIT_ASSERT(int_list->is_empty());
CPPUNIT_ASSERT_EQUAL(0, int_list->get_nb_elts());
CPPUNIT_ASSERT_EQUAL((void*)NULL, (void*)int_list->get_first());
CPPUNIT_ASSERT_EQUAL((void*)NULL, (void*)int_list->get_last());
}
void Test_ae_list::basic_tests2(void)
{
// Check get_object(int32_t pos)
CPPUNIT_ASSERT_EQUAL((void*)NULL, (void*)int_list->get_object(-1));
CPPUNIT_ASSERT_EQUAL((void*)NULL, (void*)int_list->get_object(INT_LIST_SIZE));
for (int i = 0 ; i < INT_LIST_SIZE ; i++)
{
CPPUNIT_ASSERT_EQUAL(i + 1, *(int_list->get_object(i)));
}
// Check get_node(int32_t pos)
CPPUNIT_ASSERT_EQUAL((void*)NULL, (void*)int_list->get_node(-1));
CPPUNIT_ASSERT_EQUAL((void*)NULL, (void*)int_list->get_node(INT_LIST_SIZE));
for (int i = 0 ; i < INT_LIST_SIZE-1 ; i++)
{
CPPUNIT_ASSERT_EQUAL(i + 1, *(int_list->get_node(i)->get_obj()));
}
// Check remove(node) (remove first, last and arbitrary elt)
int_list->remove(int_list->get_node(0), true, true);
int_list->remove(int_list->get_node(INT_LIST_SIZE-2), true, true);
int_list->remove(int_list->get_node(INT_LIST_SIZE/2), true, true);
ae_list<int*>* expected = new ae_list<int*>();
for (int i = 0 ; i < INT_LIST_SIZE/2 ; i++)
{
expected->add(new int(i + 2));
}
for (int i = INT_LIST_SIZE/2 ; i < INT_LIST_SIZE-3 ; i++)
{
expected->add(new int(i + 3));
}
assert_equal(expected, int_list, CPPUNIT_SOURCELINE());
expected->erase(true);
delete expected;
}
void Test_ae_list::test_extract_sublist(void)
{
// Construct the same list as int_list
ae_list<int*>* expected = new ae_list<int*>();
for (int i = 0 ; i < INT_LIST_SIZE ; i++)
{
expected->add(new int(i + 1));
}
// Initial check
assert_equal(expected, int_list, CPPUNIT_SOURCELINE());
//printf("HERE %d %d\n", expected->get_nb_elts(), int_list->get_nb_elts());
// **************************************************************************
// Extract the first element
ae_list<int*>* int_list2 = int_list->extract_sublist(0, 1);
ae_list<int*>* expected2 = new ae_list<int*>();
expected->remove(expected->get_first(), true, true);
expected2->add(new int(1));
assert_equal(expected, int_list, CPPUNIT_SOURCELINE());
assert_equal(expected2, int_list2, CPPUNIT_SOURCELINE());
expected2->erase(true);
int_list2->erase(true);
delete int_list2;
// **************************************************************************
// Extract elements 10 through 13
int_list2 = int_list->extract_sublist(10, 4);
for (int i = 0 ; i < 4 ; i++)
{
expected->remove(expected->get_node(10), true, true);
expected2->add(new int(12 + i));
}
assert_equal(expected, int_list, CPPUNIT_SOURCELINE());
assert_equal(expected2, int_list2, CPPUNIT_SOURCELINE());
expected2->erase(true);
int_list2->erase(true);
delete int_list2;
// **************************************************************************
// Extract last 4 elements
int_list2 = int_list->extract_ending_sublist(4);
for (int i = 0 ; i < 4 ; i++)
{
expected->remove(expected->get_last(), true, true);
expected2->add(new int(17 + i));
}
assert_equal(expected, int_list, CPPUNIT_SOURCELINE());
assert_equal(expected2, int_list2, CPPUNIT_SOURCELINE());
expected2->erase(true);
int_list2->erase(true);
delete int_list2;
// **************************************************************************
// Extract first 3 elements
int_list2 = int_list->extract_starting_sublist(3);
for (int i = 0 ; i < 3 ; i++)
{
expected2->add(new int(*expected->get_first()->get_obj()));
expected->remove(expected->get_first(), true, true);
}
assert_equal(expected, int_list, CPPUNIT_SOURCELINE());
assert_equal(expected2, int_list2, CPPUNIT_SOURCELINE());
expected2->erase(true);
int_list2->erase(true);
delete int_list2;
}
// ===========================================================================
// Protected Methods
// ===========================================================================
template <typename T>
void Test_ae_list::assert_equal(const ae_list<T>* expected,
const ae_list<T>* actual,
SourceLine SL)
{
// Build message string
char* msg = new char[256];
sprintf(msg, "From %s:%d", SL.fileName().c_str(), SL.lineNumber());
CPPUNIT_ASSERT_EQUAL_MESSAGE(msg,
expected->get_nb_elts(),
actual->get_nb_elts());
ae_list_node<T>* node1 = expected->get_first();
ae_list_node<T>* node2 = actual->get_first();
int32_t nb_elts = 0;
while (node1 != NULL && node2 != NULL)
{
CPPUNIT_ASSERT_EQUAL_MESSAGE(msg, *node1->get_obj(), *node2->get_obj());
nb_elts++;
node1 = node1->get_next();
node2 = node2->get_next();
}
CPPUNIT_ASSERT_EQUAL_MESSAGE(msg, expected->get_nb_elts(), nb_elts);
delete msg;
}
void Test_ae_list::testfalse() {
CPPUNIT_ASSERT(false);
}
// ===========================================================================
// Non inline accessors
// ===========================================================================
} // namespace aevol
src/libaevol/tests/run_tests.cpp deleted 100644 → 0
View file @ 5d97e536
// ****************************************************************************
//
// Aevol - An in silico experimental evolution platform
//
// ****************************************************************************
//
// Copyright: See the AUTHORS file provided with the package or <www.aevol.fr>
// Web: http://www.aevol.fr/
// E-mail: See <http://www.aevol.fr/contact/>
// Original Authors : Guillaume Beslon, Carole Knibbe, David Parsons
//
// This program is free software: you can redistribute it and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation, either version 2 of the License, or
// (at your option) any later version.
//
// This program is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
// GNU General Public License for more details.
//
// You should have received a copy of the GNU General Public License
// along with this program. If not, see <http://www.gnu.org/licenses/>.
//
//*****************************************************************************
// =================================================================
// Libraries
// =================================================================
#include <cppunit/TestCase.h>
#include <cppunit/CompilerOutputter.h>
#include <cppunit/extensions/TestFactoryRegistry.h>
#include <cppunit/ui/text/TestRunner.h>
// =================================================================
// Project Files
// =================================================================
#include "Test_JumpingMT.h"
using namespace aevol;
// ===========================================================================
// Declare Used Namespaces
// ===========================================================================
using namespace CppUnit;
// ===========================================================================
// Declare Miscellaneous Functions
// ===========================================================================
int main(int argc, char* argv[])
{
// Print message
cout << "Running regression tests";
// Get the top level suite from the registry
Test *suite = TestFactoryRegistry::getRegistry().makeTest();
// Adds the test to the list of test to run
TextUi::TestRunner runner;
runner.addTest(suite);
// Change the default outputter to a compiler error format outputter
runner.setOutputter(new CompilerOutputter(&runner.result(), cerr));
// Run the tests.
bool wasSucessful = runner.run();
// Return error code 1 if the one of test failed.
return wasSucessful ? 0 : 1;
}
// ===========================================================================
// Define Miscellaneous Functions
// ===========================================================================
src/post_treatments/CMakeLists.txt
View file @ 7ffdc5fd
......@@ -33,6 +33,9 @@ target_link_libraries(lineage PUBLIC aevol)
add_executable(mutagenesis mutagenesis.cpp)
target_link_libraries(mutagenesis PUBLIC aevol)
add_executable(mutationalrobustness mutationalrobustness.cpp)
target_link_libraries(mutationalrobustness PUBLIC aevol)
# ============================================================================
# Set behaviour on make install
......
src/post_treatments/Makefile.am
View file @ 7ffdc5fd
......@@ -39,10 +39,10 @@ aevol_misc_ancstats \
aevol_misc_create_eps \
aevol_misc_extract \
aevol_misc_lineage \
aevol_misc_mutagenesis
# aevol_misc_fixed_mutations
aevol_misc_mutagenesis \
aevol_misc_mutational_robustness \
aevol_misc_fixed_mutations
# aevol_misc_gene_families
# aevol_misc_mutational_robustness
# aevol_misc_robustness
# aevol_misc_transform_plasmid
......@@ -56,11 +56,11 @@ CLEANFILES = $(bin_PROGRAMS)
aevol_misc_ancstats_SOURCES = ancstats.cpp
aevol_misc_create_eps_SOURCES = create_eps.cpp
aevol_misc_extract_SOURCES = extract.cpp
# aevol_misc_fixed_mutations_SOURCES = fixed_mutations.cpp
aevol_misc_fixed_mutations_SOURCES = fixed_mutations.cpp
# aevol_misc_gene_families_SOURCES = gene_families.cpp
aevol_misc_lineage_SOURCES = lineage.cpp
aevol_misc_mutagenesis_SOURCES = mutagenesis.cpp
# aevol_misc_mutational_robustness_SOURCES = mutationalrobustness.cpp
aevol_misc_mutational_robustness_SOURCES = mutationalrobustness.cpp
# aevol_misc_robustness_SOURCES = robustness.cpp
# aevol_misc_template_SOURCES = template.cpp
# aevol_misc_transform_plasmid_SOURCES = transform_plasmid.cpp
......
src/post_treatments/ancstats.cpp
View file @ 7ffdc5fd
......@@ -41,7 +41,7 @@
#include <sys/stat.h>
#include <unistd.h>
#include <list>
#include <iostream>
......@@ -123,7 +123,7 @@ int main(int argc, char** argv)
static struct option long_options[] =
{
{"help", no_argument, NULL, 'h'},
{"version", no_argument, NULL, 'V' },
{"version", no_argument, NULL, 'V'},
{"verbose", no_argument, NULL, 'v'},
{"nocheck", no_argument, NULL, 'n'},
{"fullcheck", no_argument, NULL, 'c'},
......@@ -287,13 +287,12 @@ int main(int argc, char** argv)
// Prepare the initial ancestor and write its stats
// ==================================================
GridCell* grid_cell = new GridCell(lineage_file, exp_manager, nullptr);
// Individual*indiv = Individual::CreateIndividual(exp_manager, lineage_file);
auto* indiv = grid_cell->individual();
indiv->Evaluate();
indiv->compute_statistical_data();
indiv->compute_non_coding();
mystats->write_statistics_of_this_indiv(indiv);
mystats->write_statistics_of_this_indiv(indiv, nullptr);
// Optional outputs
......@@ -444,7 +443,7 @@ int main(int argc, char** argv)
indiv->compute_statistical_data();
indiv->compute_non_coding();
mystats->write_statistics_of_this_indiv(indiv);
mystats->write_statistics_of_this_indiv(indiv, rep);
// Optional outputs
write_environment_stats(time(), phenotypicTargetHandler, env_output_file);
......
src/post_treatments/fixed_mutations.cpp
View file @ 7ffdc5fd
......@@ -207,9 +207,21 @@ int main(int argc, char** argv)
// =========================
// Open the experiment manager
ae_exp_manager* exp_manager = new ae_exp_manager();
ExpManager* exp_manager = new ExpManager();
exp_manager->load(t0, true, false);
Environment* env = new Environment(*(exp_manager->env())); // independent copy
// The current version doesn't allow for phenotypic variation nor for
// different phenotypic targets among the grid
if (not exp_manager->world()->phenotypic_target_shared())
Utils::ExitWithUsrMsg("sorry, fixmut has not yet been implemented "
"for per grid-cell phenotypic target");
auto phenotypicTargetHandler =
exp_manager->world()->phenotypic_target_handler();
if (not (phenotypicTargetHandler->var_method() == NO_VAR))
Utils::ExitWithUsrMsg("sorry,fixmut has not yet been implemented "
"for variable phenotypic targets");
// GB Environment* env = new Environment(*(exp_manager->env())); // independent copy
int64_t backup_step = exp_manager->backup_step();
......@@ -263,9 +275,11 @@ int main(int argc, char** argv)
// Prepare the initial ancestor
// ==============================
ae_individual * indiv = new ae_individual(exp_manager, lineage_file);
indiv->evaluate(habitat);
GridCell* grid_cell = new GridCell(lineage_file, exp_manager, nullptr);
auto* indiv = grid_cell->individual();
indiv->Evaluate();
indiv->compute_statistical_data();
indiv->compute_non_coding();
if (verbose)
{
......@@ -281,29 +295,31 @@ int main(int argc, char** argv)
// time a backup is available)
// ===============================================================================
ae_individual* stored_indiv = NULL;
std::list<GeneticUnit>::const_iterator stored_unit;
Individual* stored_indiv = NULL;
// GB std::list<GeneticUnit>::const_iterator stored_unit;
ReplicationReport* rep = nullptr;
int32_t index, genetic_unit_number, unitlen_before;
int32_t nb_genes_at_breakpoints, nb_genes_in_segment, nb_genes_in_replaced_segment;
double metabolic_error_before, metabolic_error_after, impact_on_metabolic_error;
char mut_descr_string[80];
ae_exp_manager* exp_manager_backup = NULL;
Environment* backup_env = NULL;
ExpManager* exp_manager_backup = NULL;
Habitat* backup_habitat = NULL;
bool check_now = false;
aevol::Time::plusplus();
aevol::AeTime::plusplus();
while (time() <= t_end)
{
ae_replication_report* rep = new ae_replication_report(lineage_file, indiv);
rep = new ReplicationReport(lineage_file, indiv);
index = rep->id(); // who are we building...
indiv->set_replication_report(rep);
//GB indiv->set_replication_report(rep);
// Check now?
check_now = ((check == FULL_CHECK && ae_utils::mod(time(), backup_step) == 0) ||
(check == ENV_CHECK && ae_utils::mod(time(), backup_step) == 0) ||
check_now = ((check == FULL_CHECK && Utils::mod(time(), backup_step) == 0) ||
(check == ENV_CHECK && Utils::mod(time(), backup_step) == 0) ||
(check == LIGHT_CHECK && time() == t_end));
......@@ -311,104 +327,102 @@ int main(int argc, char** argv)
printf("Rebuilding ancestor at generation %" PRId64 " (index %" PRId32 ")...",
time(), index);
env->build();
/* GB env->build();
env->apply_variation();
indiv->reevaluate(env);
indiv->reevaluate(env); */
indiv->Reevaluate();
// TODO <david.parsons@inria.fr> Check for phenotypic variation has to be
// done for all the grid cells, disable checking until coded
// Check, and possibly update, the environment according to the backup files (update necessary if the env. was modified by aevol_modify at some point)
if (ae_utils::mod(time(), backup_step) == 0)
{
char env_file_name[255];
sprintf(env_file_name, "./" ENV_FNAME_FORMAT, time());
gzFile env_file = gzopen(env_file_name, "r");
backup_env = new Environment();
backup_env->load(env_file);
if (! env->is_identical_to(*backup_env, tolerance))
{
printf("Warning: At t=%" PRId64 ", the replayed environment is not the same\n", time());
printf(" as the one saved at generation %" PRId64 "... \n", time());
printf(" with tolerance of %lg\n", tolerance);
printf("Replacing the replayed environment by the one stored in the backup.\n");
delete env;
env = new Environment(*backup_env);
}
delete backup_env;
}
// if (Utils::mod(time(), backup_step) == 0)
// {
// char env_file_name[255];
// sprintf(env_file_name, "./" ENV_FNAME_FORMAT, time());
// gzFile env_file = gzopen(env_file_name, "r");
// backup_habitat = new Habitat();
// backup_habitat->load(env_file);
//
// if (! env->is_identical_to(*backup_habitat, tolerance))
// {
// printf("Warning: At t=%" PRId64 ", the replayed environment is not the same\n", time());
// printf(" as the one saved at generation %" PRId64 "... \n", time());
// printf(" with tolerance of %lg\n", tolerance);
// printf("Replacing the replayed environment by the one stored in the backup.\n");
// delete env;
// env = new Habitat(*backup_habitat);
// }
// delete backup_habitat;
// }
// Warning: this portion of code won't work if the number of units changes
// during the evolution, or if some translocations occurred between different genetic units
genetic_unit_number = 0;
std::list<DnaReplicReport*>::const_iterator dnareport = rep->dna_replic_reports().begin();
std::list<GeneticUnit>::iterator unit = indiv->genetic_unit_list_nonconst().begin();
// GB genetic_unit_number = 0;
// std::list<DnaReplicationReport*>::const_iterator dnareport = rep->dna_replic_report().begin();
// GB std::list<GeneticUnit>::iterator unit = indiv->genetic_unit_list_nonconst().begin(); */
if (check_now && ae_utils::mod(time(), backup_step) == 0)
GeneticUnit& gen_unit = indiv->genetic_unit_nonconst(0);
GeneticUnit* stored_gen_unit = nullptr;
Individual* stored_indiv = nullptr;
if (check_now && Utils::mod(time(), backup_step) == 0)
{
exp_manager_backup = new ae_exp_manager();
exp_manager_backup = new ExpManager();
exp_manager_backup->load(time(), true, false);
// TODO: disabled tmp
// stored_indiv = new ae_individual(* (ae_individual *)exp_manager_backup->indiv_by_id(index), false);
stored_unit = stored_indiv->genetic_unit_list().begin();
stored_indiv = new Individual(
*(Individual*) exp_manager_backup->indiv_by_id(index));
stored_gen_unit = &(stored_indiv->genetic_unit_nonconst(0));
}
while (dnareport != rep->dna_replic_reports().end())
{
assert(unit != indiv->genetic_unit_list().end());
unit->dna()->set_replic_report(*dnareport);
// ***************************************
// Transfer events
// ***************************************
for (const auto& mutation: (*dnareport)->HT()) {
metabolic_error_before = indiv->dist_to_target_by_feature(METABOLISM);
unitlen_before = unit->dna()->length();
unit->compute_nb_of_affected_genes(&mutation, nb_genes_at_breakpoints, nb_genes_in_segment, nb_genes_in_replaced_segment);
unit->dna()->undergo_this_mutation(&mutation);
indiv->reevaluate(env);
metabolic_error_after = indiv->dist_to_target_by_feature(METABOLISM);
impact_on_metabolic_error = metabolic_error_after - metabolic_error_before;
mutation.generic_description_string(mut_descr_string);
fprintf(output,
"%" PRId64 " %" PRId32 " %s %" PRId32 " %.15f %" PRId32
" %" PRId32 " %" PRId32 " \n",
time(), genetic_unit_number, mut_descr_string, unitlen_before,
impact_on_metabolic_error, nb_genes_at_breakpoints, nb_genes_in_segment, nb_genes_in_replaced_segment);
}
// For each genetic unit, replay the replication (undergo all mutations)
// TODO <david.parsons@inria.fr> disabled for multiple GUs
const auto& dnarep = rep->dna_replic_report();
// TODO <guillaume.beslon@inria.fr> compute statistics for HT
// for (const auto& mut: dnarep.HT()) {
// metabolic_error_before = indiv->dist_to_target_by_feature(METABOLISM);
// unitlen_before = gen_unit.dna()->length();
// gen_unit.compute_nb_of_affected_genes(*mut, nb_genes_at_breakpoints, nb_genes_in_segment, nb_genes_in_replaced_segment);
// gen_unit.dna()->undergo_this_mutation(*mut);
// indiv->Reevaluate(env);
// metabolic_error_after = indiv->dist_to_target_by_feature(METABOLISM);
// impact_on_metabolic_error = metabolic_error_after - metabolic_error_before;
// mutation.generic_description_string(mut_descr_string);
// fprintf(output,
// "%" PRId64 " %" PRId32 " %s %" PRId32 " %.15f %" PRId32
// " %" PRId32 " %" PRId32 " \n",
// time(), genetic_unit_number, mut_descr_string, unitlen_before,
// impact_on_metabolic_error, nb_genes_at_breakpoints, nb_genes_in_segment, nb_genes_in_replaced_segment);
// }
// ***************************************
// Rearrangement events
// ***************************************
for (const auto& mutation: (*dnareport)->rearrangements()) {
for (const auto& mutation: dnarep.rearrangements()) {
metabolic_error_before = indiv->dist_to_target_by_feature(METABOLISM);
unitlen_before = unit->dna()->length();
unit->compute_nb_of_affected_genes(&mutation, nb_genes_at_breakpoints, nb_genes_in_segment, nb_genes_in_replaced_segment);
unitlen_before = gen_unit.dna()->length();
//GB gen_unit.compute_nb_of_affected_genes(*mutation, nb_genes_at_breakpoints, nb_genes_in_segment, nb_genes_in_replaced_segment);
unit->dna()->undergo_this_mutation(&mutation);
gen_unit.dna()->undergo_this_mutation(*mutation);
indiv->reevaluate(env);
indiv->Reevaluate();
metabolic_error_after = indiv->dist_to_target_by_feature(METABOLISM);
impact_on_metabolic_error = metabolic_error_after - metabolic_error_before;
mutation.generic_description_string(mut_descr_string);
mutation->generic_description_string(mut_descr_string);
fprintf(output,
"%" PRId64 " %" PRId32 " %s %" PRId32 " %.15f %" PRId32
" %" PRId32 " %" PRId32 " \n",
time(), genetic_unit_number, mut_descr_string, unitlen_before,
impact_on_metabolic_error, nb_genes_at_breakpoints, nb_genes_in_segment,
nb_genes_in_replaced_segment);
"%" PRId64 " %" PRId32 " %s %" PRId32 " %.15f \n",
time(), 0, mut_descr_string, unitlen_before,
impact_on_metabolic_error);
}
......@@ -416,100 +430,85 @@ int main(int argc, char** argv)
// Local events (point mutations & small indels)
// ***************************************
for (const auto& mutation: (*dnareport)->mutations()) {
for (const auto& mutation: dnarep.mutations()) {
metabolic_error_before = indiv->dist_to_target_by_feature(METABOLISM);
unitlen_before = unit->dna()->length();
unit->compute_nb_of_affected_genes(&mutation, nb_genes_at_breakpoints, nb_genes_in_segment, nb_genes_in_replaced_segment);
unitlen_before = gen_unit.dna()->length();
// GB gen_unit.compute_nb_of_affected_genes(&mutation, nb_genes_at_breakpoints, nb_genes_in_segment, nb_genes_in_replaced_segment);
unit->dna()->undergo_this_mutation(&mutation);
gen_unit.dna()->undergo_this_mutation(*mutation);
indiv->reevaluate(env);
indiv->Reevaluate();
metabolic_error_after = indiv->dist_to_target_by_feature(METABOLISM);
impact_on_metabolic_error = metabolic_error_after - metabolic_error_before;
mutation.generic_description_string(mut_descr_string);
mutation->generic_description_string(mut_descr_string);
fprintf(output,
"%" PRId64 " %" PRId32 " %s %" PRId32 " %.15f %" PRId32
" %" PRId32 " %" PRId32 " \n",
time(), genetic_unit_number, mut_descr_string, unitlen_before,
impact_on_metabolic_error, nb_genes_at_breakpoints, nb_genes_in_segment,
nb_genes_in_replaced_segment);
"%" PRId64 " %" PRId32 " %s %" PRId32 " %.15f \n",
time(), 0, mut_descr_string, unitlen_before,
impact_on_metabolic_error);
}
if (check_now && ae_utils::mod(time(), backup_step) == 0)
if (check_now)
{
if (verbose)
{
if (verbose)
{
printf("Checking the sequence of the unit...");
fflush(NULL);
}
assert(stored_unit != stored_indiv->genetic_unit_list().end());
printf("Checking the sequence of the unit...");
fflush(NULL);
}
char * str1 = new char[unit->dna()->length() + 1];
memcpy(str1, unit->dna()->data(), \
unit->dna()->length()*sizeof(char));
str1[unit->dna()->length()] = '\0';
char * str1 = new char[gen_unit.dna()->length() + 1];
memcpy(str1, gen_unit.dna()->data(), \
gen_unit.dna()->length()*sizeof(char));
str1[gen_unit.dna()->length()] = '\0';
char * str2 = new char[(stored_unit->dna())->length() + 1];
memcpy(str2, (stored_unit->dna())->data(), (stored_unit->dna())->length()*sizeof(char));
str2[(stored_unit->dna())->length()] = '\0';
char * str2 = new char[(stored_gen_unit->dna())->length() + 1];
memcpy(str2, (stored_gen_unit->dna())->data(),
(stored_gen_unit->dna())->length()*sizeof(char));
str2[(stored_gen_unit->dna())->length()] = '\0';
if(strncmp(str1,str2, (stored_unit->dna())->length())==0)
{
if (verbose) printf(" OK\n");
}
else
{
if (verbose) printf(" ERROR !\n");
fprintf(stderr, "Error: the rebuilt unit is not the same as \n");
fprintf(stderr, "the one saved at generation %" PRId64 "... ", time());
fprintf(stderr, "Rebuilt unit : %zu bp\n %s\n", strlen(str1), str1);
fprintf(stderr, "Stored unit : %zu bp\n %s\n", strlen(str2), str2);
delete [] str1;
delete [] str2;
gzclose(lineage_file);
delete indiv;
delete stored_indiv;
delete exp_manager_backup;
delete exp_manager;
exit(EXIT_FAILURE);
}
if (strncmp(str1, str2, stored_gen_unit->dna()->length()) == 0) {
if (verbose)
printf(" OK\n");
}
else {
if (verbose) printf(" ERROR !\n");
fprintf(stderr, "Error: the rebuilt genetic unit is not the same as \n");
fprintf(stderr, "the one saved at generation %" PRId64 "... ", time());
fprintf(stderr, "Rebuilt unit : %" PRId32 " bp\n %s\n", (int32_t)strlen(str1), str1);
fprintf(stderr, "Stored unit : %" PRId32 " bp\n %s\n", (int32_t)strlen(str2), str2);
delete [] str1;
delete [] str2;
++stored_unit;
gzclose(lineage_file);
delete indiv;
delete stored_indiv;
delete exp_manager_backup;
delete exp_manager;
exit(EXIT_FAILURE);
}
++dnareport;
++unit;
genetic_unit_number++;
}
assert(unit == indiv->genetic_unit_list().end());
delete [] str1;
delete [] str2;
}
if (verbose) printf(" OK\n");
delete rep;
if (check_now && ae_utils::mod(time(), backup_step) == 0)
{
assert(stored_unit == stored_indiv->genetic_unit_list().end());
delete stored_indiv;
delete exp_manager_backup;
}
aevol::Time::plusplus();
if (check_now)
{
delete stored_indiv;
delete exp_manager_backup;
}
aevol::AeTime::plusplus();
}
gzclose(lineage_file);
fclose(output);
delete exp_manager;
delete indiv;
delete env;
exit(EXIT_SUCCESS);
......
src/post_treatments/lineage.cpp
View file @ 7ffdc5fd
......@@ -368,6 +368,7 @@ int main(int argc, char** argv)
// NB : The list of individuals is sorted according to the index
const Individual* const initial_ancestor = exp_manager->indiv_by_id(indices[0]);
// Write file "header"
gzwrite(lineage_file, &t0, sizeof(t0));
gzwrite(lineage_file, &t_end, sizeof(t_end));
gzwrite(lineage_file, &final_indiv_index, sizeof(final_indiv_index));
......
src/post_treatments/mutationalrobustness.cpp
View file @ 7ffdc5fd
......@@ -32,140 +32,314 @@
of a population.
*/
#include <list>
#include <getopt.h>
#include <libgen.h>
#include <list>
#include "aevol.h"
using namespace aevol;
enum check_type {
FULL_CHECK = 0,
LIGHT_CHECK = 1,
ENV_CHECK = 2,
NO_CHECK = 3
};
void print_help(char* prog_name);
void analyse_indiv(ae_exp_manager*, ae_individual*, FILE*, int32_t);
void analyze_indiv(Individual* initial_indiv,
FILE* output,
int32_t ndiv,
std::shared_ptr<JumpingMT> prng,
bool verbose);
int main(int argc, char* argv[])
{
int main(int argc, char* argv[]) {
// Load parameters from command line
int32_t ndiv = 100000; // Default number of mutants per individual
int32_t gener = -1; // What generation to load
char* output_file_name = NULL;
bool best_only = false; // Treat only the best individual?
// TODO <david.parsons@inria.fr> version
const char * options_list = "hn:r:o:b";
static struct option long_options_list[] = {
{ "help", 0, NULL, 'h' },
{ "number", 1, NULL, 'n' },
{ "generation", 1, NULL, 'r' },
{ "output", 1, NULL, 'o' },
{ "best", 0, NULL, 'b' },
{ 0, 0, 0, 0 }
int32_t nb_indiv = 1000; // Default number of mutants per individual
int32_t begin = 0; // Default starting generation
int32_t end = -1; // Default ending generation (-1 for last generation stored in lineage file)
int32_t period = 1; // Period of analyze
char* output_file_name = NULL;
char* lineage_file_name = NULL;
bool verbose = false;
const char* short_options = "hVvn:b:e:p:f:o:l";
static struct option long_options[] =
{
{"help", no_argument, NULL, 'h'},
{"version", no_argument, NULL, 'V'},
{"verbose", no_argument, NULL, 'v'},
{"number", required_argument, NULL, 'n'},
{"begin", required_argument, NULL, 'b'},
{"end", required_argument, NULL, 'e'},
{"period", required_argument, NULL, 'p'},
{"file", required_argument, NULL, 'f'},
{"output", required_argument, NULL, 'o'},
{0, 0, 0, 0}
};
int option;
while ((option = getopt_long(argc, argv, options_list, long_options_list, NULL)) != -1)
{
switch (option)
{
case 'h' : print_help(basename(argv[0])); exit(EXIT_SUCCESS); break;
case 'n' : ndiv = atol(optarg); break;
case 'r' : gener = atol(optarg); break;
while ((option = getopt_long(argc, argv, short_options, long_options, NULL)) != -1) {
switch (option) {
case 'h' :
print_help(argv[0]);
exit(EXIT_SUCCESS);
case 'V' :
Utils::PrintAevolVersion();
exit(EXIT_SUCCESS);
case 'v' :
verbose = true;
break;
case 'n' :
nb_indiv = atol(optarg);
break;
case 'b' :
begin = atol(optarg);
break;
case 'e' :
end = atol(optarg);
break;
case 'p' :
period = atol(optarg);
break;
case 'o' :
{
output_file_name = new char[strlen(optarg) + 1];
sprintf(output_file_name, "%s", optarg);
break;
}
case 'b' : best_only = true; break;
case 'f' :
lineage_file_name = new char[strlen(optarg) + 1];
sprintf(lineage_file_name, "%s", optarg);
break;
}
}
// Load the population from the backup file
if (gener == -1){
printf("You must specify a generation number. Please use the option -r or --generation.\n");
// =======================
// Open the lineage file
// =======================
gzFile lineage_file = gzopen(lineage_file_name, "r");
if (lineage_file == Z_NULL) {
fprintf(stderr, "ERROR : Could not read lineage file %s\n", lineage_file_name);
exit(EXIT_FAILURE);
}
ae_exp_manager* exp_manager = new ae_exp_manager();
exp_manager->load(gener, true, false);
// int32_t nb_indivs = exp_manager->nb_indivs();
int64_t t0 = 0;
int64_t t_end = 0;
int32_t final_indiv_index = 0;
int32_t final_indiv_rank = 0;
// Open output file and write the header
FILE * output = fopen(output_file_name, "w");
if (output == NULL){
fprintf(stderr, "ERROR : Could not create the output file %s\n", output_file_name);
exit(EXIT_FAILURE);
}
gzread(lineage_file, &t0, sizeof(t0));
gzread(lineage_file, &t_end, sizeof(t_end));
gzread(lineage_file, &final_indiv_index, sizeof(final_indiv_index));
gzread(lineage_file, &final_indiv_rank, sizeof(final_indiv_rank));
// Positive impact means
fprintf(output, "# #################################################################\n");
fprintf(output, "# Mutations produced by mutationalrobustness\n");
fprintf(output, "# #################################################################\n");
fprintf(output, "# Number of replicate per individual : %" PRId32 " \n",ndiv);
fprintf(output, "# Impact on metabolism SPACE impact on secretion\n");
fprintf(output, "#\n");
// Parse and treat the individuals
if (!best_only){
for (ae_individual* indiv: exp_manager->world()->indivs())
analyse_indiv(exp_manager, indiv, output, ndiv);
if (verbose) {
printf("\n\n");
printf("===============================================================================\n");
printf(" Robustness of the ancestors of indiv. %" PRId32
" (rank %" PRId32 ") from time %" PRId64 " to %" PRId64 "\n",
final_indiv_index, final_indiv_rank, t0, t_end);
printf("================================================================================\n");
}
else{
ae_individual* indiv = exp_manager->world()->best_indiv();
analyse_indiv(exp_manager, indiv, output, ndiv);
// =============================
// Open the experience manager
// =============================
ExpManager* exp_manager = new ExpManager();
exp_manager->load(t0, true, false);
// The current version doesn't allow for phenotypic variation nor for
// different phenotypic targets among the grid
if (not exp_manager->world()->phenotypic_target_shared())
Utils::ExitWithUsrMsg("sorry, ancestor stats has not yet been implemented "
"for per grid-cell phenotypic target\n");
auto phenotypicTargetHandler =
exp_manager->world()->phenotypic_target_handler();
if (phenotypicTargetHandler->var_method() != NO_VAR)
Utils::ExitWithUsrMsg("sorry, ancestor stats has not yet been implemented "
"for variable phenotypic targets\n");
int64_t backup_step = exp_manager->backup_step();
// =========================
// Open the output file(s)
// =========================
// // Create missing directories
// int status;
// status = mkdir("stats/ancstats/", 0755);
// if ((status == -1) && (errno != EEXIST))
// err(EXIT_FAILURE, "stats/ancstats/");
FILE * output_summary = fopen(output_file_name, "w");
std::shared_ptr<JumpingMT> prng = std::make_shared<JumpingMT>(9695);
// ==============================
// Prepare the initial ancestor
// ==============================
GridCell* grid_cell = new GridCell(lineage_file, exp_manager, nullptr);
auto* indiv = grid_cell->individual();
indiv->Evaluate();
// indiv->compute_statistical_data();
// indiv->compute_non_coding();
// ==============================
// Compute robustness of the initial ancestor
// ==============================
if (begin == 0)
analyze_indiv(indiv,output_summary,nb_indiv,prng,verbose);
// ==========================================================================
// Replay the mutations to get the successive ancestors and analyze them
// ==========================================================================
ReplicationReport* rep = nullptr;
int32_t index;
ExpManager* exp_manager_backup = nullptr;
Habitat *backup_habitat = nullptr;
aevol::AeTime::plusplus();
while ((time() <= t_end) && (((time() < end)||(end == -1))))
{
rep = new ReplicationReport(lineage_file, indiv);
index = rep->id(); // who we are building...
indiv->Reevaluate();
if (verbose)
printf("Ancestor at generation %" PRId64
" has index %" PRId32 "\n", time(), index);
// 2) Replay replication (create current individual's child)
GeneticUnit& gen_unit = indiv->genetic_unit_nonconst(0);
GeneticUnit* stored_gen_unit = nullptr;
Individual* stored_indiv = nullptr;
// For each genetic unit, replay the replication (undergo all mutations)
// TODO <david.parsons@inria.fr> disabled for multiple GUs
const auto& dnarep = rep->dna_replic_report();
for (const auto& mut: dnarep.HT())
gen_unit.dna()->undergo_this_mutation(*mut);
for (const auto& mut: dnarep.rearrangements())
gen_unit.dna()->undergo_this_mutation(*mut);
for (const auto& mut: dnarep.mutations())
gen_unit.dna()->undergo_this_mutation(*mut);
// 3) All the mutations have been replayed, we can now evaluate the new individual
indiv->Reevaluate();
// if we are between "begin" and "end" and at the correct period, compute robustness
if ((time() >= begin)&&((time() < end)||(end == -1)))
{
if (((time() - begin) % period) == 0)
{
analyze_indiv(indiv,output_summary,nb_indiv,prng,verbose);
}
}
delete rep;
aevol::AeTime::plusplus();
}
// Clean memory and exit
gzclose(lineage_file);
fclose(output_summary);
delete exp_manager;
delete [] output_file_name;
// delete mystats;
delete indiv;
return EXIT_SUCCESS;
}
// Treatment of one individual
void analyse_indiv(ae_exp_manager* exp, ae_individual* initial_indiv, FILE* output, int32_t ndiv){
Environment* env = exp->env();
double initial_metabolic_error = initial_indiv->dist_to_target_by_feature(METABOLISM);
double initial_secretion_error = initial_indiv->dist_to_target_by_feature(SECRETION);
double final_metabolic_error = 0.0;
double impact_on_metabolic_error = 0.0;
double final_secretion_error = 0.0;
double impact_on_secretion_error = 0.0;
ae_individual* indiv = NULL;
int32_t i;
// Perform ndiv reproductions with mutations
for (i = 0; i < ndiv; i++){
if (i % 1000 == 0){
printf("*");
fflush(stdout);
void analyze_indiv(Individual* indiv,
FILE* output_summary,
int32_t nb_indiv,
std::shared_ptr<JumpingMT> prng,
bool verbose)
{
double fracplus = 0;
double meanplus = 0;
double fracminus = 0;
double meanminus = 0;
double maxplus = 0;
double maxminus = 0;
double fracnull = 0;
int32_t nb_events = 0;
double indiv_metabolic_error = indiv->dist_to_target_by_feature(METABOLISM); const Habitat& habitat = indiv->habitat();
for (int32_t i = 0; i < nb_indiv; i++)
{
Individual* new_indiv = new Individual(indiv, 0, prng, prng);
// Perform transfer, rearrangements and mutations
if (not new_indiv->allow_plasmids()) {
const GeneticUnit* chromosome = &new_indiv->genetic_unit_list().front();
nb_events = chromosome->dna()->perform_mutations(indiv->id());
}
else {
printf("WARNING: Mutational Robustness does not handle multiple Genetic Units\n");
}
if (nb_events == 0)
{
fracnull++;
}
else
{
new_indiv->EvaluateInContext(habitat);
double new_metabolic_error = new_indiv->dist_to_target_by_feature(METABOLISM);
if (new_metabolic_error == indiv_metabolic_error)
fracnull++;
if (new_metabolic_error > indiv_metabolic_error)
{
fracminus++;
if ((new_metabolic_error - indiv_metabolic_error) > maxminus) maxminus = new_metabolic_error - indiv_metabolic_error;
meanminus += new_metabolic_error - indiv_metabolic_error;
}
if (new_metabolic_error < indiv_metabolic_error)
{
fracplus++;
if ((new_metabolic_error - indiv_metabolic_error) < maxplus) maxplus = new_metabolic_error - indiv_metabolic_error;
meanplus += new_metabolic_error - indiv_metabolic_error;
}
}
delete new_indiv;
}
indiv = exp->sel()->do_replication(initial_indiv, -1);
indiv->reevaluate(env);
final_metabolic_error = indiv->dist_to_target_by_feature(METABOLISM);
impact_on_metabolic_error = final_metabolic_error - initial_metabolic_error;
final_secretion_error = indiv->dist_to_target_by_feature(SECRETION);
impact_on_secretion_error = final_secretion_error - initial_secretion_error;
fprintf(output, "%+.15f %+.15f \n",impact_on_metabolic_error, impact_on_secretion_error);
delete indiv;
if (verbose) printf("f+: %f f0: %f f-:%f\n",fracplus,fracnull,fracminus);
if ( output_summary == NULL ){
fprintf( stderr, "ERROR : Could not create robustness_summary.txt\n");
exit( EXIT_FAILURE );
}
fprintf(output, "\n");
fprintf(output_summary, "%lld %.15f %.15f %.15f %.15f %.15f %.15f %.15f\n",time(),fracplus/nb_indiv,fracnull/nb_indiv,fracminus/nb_indiv, meanplus/fracplus,meanminus/fracminus,maxplus,maxminus);
}
// Print help
void print_help(char* prog_name){
void print_help(char* prog_name) {
printf("\n\
%s is a post-treatment that generates and analyses a large quantity of mutants for all individuals in a backup. For each mutants we record the phenotypic effect on metabolism and on secretion.\n\n\
Usage: %s [-h] -r num_generation -o output_file_name -n num_mutants [-b] \n\
\t-h : display this screen\n\
\t-r num_generation : read the generation num_generation from a full aevol backup\n\
\t-o output_file_name : write the results in file output_file_name\n\
\t-n num_mutants : generate and analyse num_mutants per individual\n\
\t-b : only treat the best individual\n\n\
Example:\n\t%s -r 20000 -n 1000 -o toto.out\n",prog_name,prog_name,prog_name);
%s is a post-treatment that generates and analyses a large quantity of mutants for a lineage of ancestors.\
For each mutant we record the phenotypic effect on metabolism.\n\n\
Usage: %s [-h] -i input_file_name -o output_file_name [-b start_at_generation] [-e end_at_generation] [-p period] [-n num_mutants] [-r] [-h bin_size] [-v verbose] [-s mutation_seed]\n\
\t-h: display this screen\n\
\t-f input_file_name: lineage file to be analyzed\n\
\t-o output_file_name: name of the output file (to be written in ./stats/ancstats). In case of histogram output (-h) one file will be produced for each histogram and output_file_name will be postfixed with the generation number\n\
\t-b start_at_generation: first generation of the lineage to be analyzed (default: 0)\n\
\t-e end_at_generation: last generation of the lineage to be analyzed (default: last generation stored in the input file)\n\
\t-p period: temporal resolution of the analyze (default: 1)\n\
\t-n nb_mutants : generate and analyse nb_mutants per individual (default: 1000)\n\
Example:\n\t%s -i lineage_file -o toto.out -b 4000 -e 5000 -p 10 -n 100000 -s 19769\n", prog_name, prog_name, prog_name);
// \t-r: raw output; store the difference of metabolic error for each mutant generated (warning: the output file may quickly grow)\n
// \t-h bin_size: store the histogram with a bin_size resolution. One output file is generated for each histogram (postfixed with the generation number)\n
// \t-s mutation_seed: specify the seed to be used for the mutation random generator\n\n
}
src/post_treatments/robustness.cpp
View file @ 7ffdc5fd
......@@ -357,10 +357,10 @@ void print_help(char* prog_path)
printf(" or : %s -g GENER [-n NBCHILDREN] [-r RANK | -i INDEX]\n", prog_name);
printf("\n");
printf("This program computes the replication statistics of all the individuals of a given generation,\n");
printf("like the proportion of neutral, beneficial, deleterious offsprings. This is done by simulating\n");
printf("like the proportion of neutral, beneficial and deleterious offsprings. This is done by simulating\n");
printf("NBCHILDREN replications for each individual, with its mutation, rearrangement and transfer rates.\n");
printf("Depending on those rates and genome size, there can be several events per replication.\n");
printf("Those statistics are written in analysis-generationGENER/robustness-allindivs-gGENER.out, with one \n");
printf("These statistics are written in analysis-generationGENER/robustness-allindivs-gGENER.out, with one \n");
printf("line per individual in the specified generation. \n\n");
printf("The program also outputs detailed statistics for one of the individuals (the best one by default). \n");
printf("The detailed statistics for this individual are written in the file called \n");
......
src/tests/.gitignore
View file @ 7ffdc5fd
# Ignore test binaries
individualTest
sampleTest
IndividualTest
TranscriptionTranslationTest
sampleTest
# ignore test log and result files
**/*.log
**/*.trs
src/tests/IndividualTest.cpp
View file @ 7ffdc5fd
......@@ -59,8 +59,7 @@ using std::list;
// Class IndividualTest #
// #
//############################################################################
class IndividualTest : public testing::Test
{
class IndividualTest : public testing::Test {
protected:
virtual void SetUp(void);
virtual void TearDown(void);
......@@ -76,123 +75,143 @@ class IndividualTest : public testing::Test
// ===========================================================================
void IndividualTest::SetUp(void)
{
// Build ad-hoc genomes
// (and reverse to test the same things on the lagging strand.):
//
// indiv1: (AS + prom + AS + AG + AS + term + AS + prom + AS)
// indiv2: reverse
// indiv3: (AS + AG + AS + term + AS + prom + AS)
// indiv4: reverse
//
// AS = Arbitrary Sequence
// AG = Arbitrary Gene
// Do not modify the sequences !
ExpSetup* exp_s = new ExpSetup(nullptr);
for (int i = 0; i <= 1; i++) {
exp_s->set_fuzzy_flavor(i);
// Define a few arbitrary sequences
char as[5][10] = {
"0011",
"11101",
"110011",
"11000",
"000101"
};
if (FuzzyFactory::fuzzyFactory == NULL)
FuzzyFactory::fuzzyFactory = new FuzzyFactory(exp_s);
else {
delete FuzzyFactory::fuzzyFactory;
FuzzyFactory::fuzzyFactory = new FuzzyFactory(exp_s);
}
printf("F %d\n", i);
// Build ad-hoc genomes
// (and reverse to test the same things on the lagging strand.):
//
// indiv1: (AS + prom + AS + AG + AS + term + AS + prom + AS)
// indiv2: reverse
// indiv3: (AS + AG + AS + term + AS + prom + AS)
// indiv4: reverse
//
// AS = Arbitrary Sequence
// AG = Arbitrary Gene
// Do not modify the sequences !
// Define an arbitrary gene
char gene[255];
sprintf(gene, "%s0011000100110110010001", SHINE_DAL_SEQ);
// Define a few arbitrary sequences
char as[5][10] = {
"0011",
"11101",
"110011",
"11000",
"000101"
};
// Define an arbitrary terminator
char term[TERM_SIZE+1] = "01000001101";
// Define an arbitrary gene
char gene[255];
sprintf(gene, "%s0011000100110110010001", SHINE_DAL_SEQ);
// Define a few arbitrary promoters
char prom[2][23] = {
"0101010001110110010110", // dist from consensus: 2 => basal level: 0.6
"0101011001110010010010" // dist from consensus: 1 => basal level: 0.8
};
// Define an arbitrary terminator
char term[TERM_SIZE + 1] = "01000001101";
// Construct a gene with these arbitrary sequences
char* genome = new char[1024];
sprintf(genome, "%s%s%s%s%s%s%s%s%s", as[0], prom[0], as[1], gene, as[2],
term, as[3], prom[1], as[4]);
// Define a few arbitrary promoters
char prom[2][23] = {
"0101010001110110010110", // dist from consensus: 2 => basal level: 0.6
"0101011001110010010010" // dist from consensus: 1 => basal level: 0.8
};
// Build indiv1
MutationParams params_mut;
indiv1 = new Individual(nullptr, nullptr, nullptr, std::make_shared<MutationParams>(params_mut), 1.0, 10, 1000, false, 1, "anon-strain-1", 0);
indiv1->add_GU(genome, strlen(genome));
genome = NULL;
// Construct a gene with these arbitrary sequences
char* genome = new char[1024];
sprintf(genome, "%s%s%s%s%s%s%s%s%s", as[0], prom[0], as[1], gene, as[2],
term, as[3], prom[1], as[4]);
// Do transcription and translation
indiv1->do_transcription();
indiv1->do_translation();
// Build indiv1
MutationParams params_mut;
indiv1 = new Individual(nullptr, nullptr, nullptr,
std::make_shared<MutationParams>(params_mut), 1.0,
10, 1000, false, 1, "anon-strain-1", 0);
indiv1->add_GU(genome, strlen(genome));
genome = NULL;
// Do transcription and translation
indiv1->do_transcription();
indiv1->do_translation();
// Build indiv2
// Reverse the whole genome
genome = indiv1->get_genetic_unit(0).get_dna()->get_subsequence(0,0,LAGGING);
indiv2 = new Individual(nullptr, nullptr, nullptr, std::make_shared<MutationParams>(params_mut), 1.0, 10, 1000, false, 1, "anon-strain-2", 0);
indiv2->add_GU(genome, strlen(genome));
genome = NULL;
// Do transcription and translation
indiv2->do_transcription();
indiv2->do_translation();
// Build indiv2
// Reverse the whole genome
genome = indiv1->genetic_unit(0).dna()->subsequence(0, 0, LAGGING);
indiv2 = new Individual(nullptr, nullptr, nullptr,
std::make_shared<MutationParams>(params_mut), 1.0,
10, 1000, false, 1, "anon-strain-2", 0);
indiv2->add_GU(genome, strlen(genome));
genome = NULL;
// Do transcription and translation
indiv2->do_transcription();
indiv2->do_translation();
// Build indiv3
genome = new char[1024];
sprintf(genome, "%s%s%s%s%s%s%s", as[0], gene, as[1], term, as[2], prom[1], as[3]);
indiv3 = new Individual(nullptr, nullptr, nullptr, std::make_shared<MutationParams>(params_mut), 1.0, 10, 1000, false, 1, "anon-strain-3", 0);
indiv3->add_GU(genome, strlen(genome));
genome = NULL;
// Do transcription and translation
indiv3->do_transcription();
indiv3->do_translation();
// Build indiv3
genome = new char[1024];
sprintf(genome, "%s%s%s%s%s%s%s", as[0], gene, as[1], term, as[2], prom[1],
as[3]);
indiv3 = new Individual(nullptr, nullptr, nullptr,
std::make_shared<MutationParams>(params_mut), 1.0,
10, 1000, false, 1, "anon-strain-3", 0);
indiv3->add_GU(genome, strlen(genome));
genome = NULL;
// Do transcription and translation
indiv3->do_transcription();
indiv3->do_translation();
// Build indiv4
genome = indiv3->get_genetic_unit(0).get_dna()->get_subsequence(0,0,LAGGING);
indiv4 = new Individual(nullptr, nullptr, nullptr, std::make_shared<MutationParams>(params_mut), 1.0, 10, 1000, false, 1, "anon-strain-4", 0);
indiv4->add_GU(genome, strlen(genome));
genome = NULL;
// Do transcription and translation
indiv4->do_transcription();
indiv4->do_translation();
// Build indiv4
genome = indiv3->genetic_unit(0).dna()->subsequence(0, 0, LAGGING);
indiv4 = new Individual(nullptr, nullptr, nullptr,
std::make_shared<MutationParams>(params_mut), 1.0,
10, 1000, false, 1, "anon-strain-4", 0);
indiv4->add_GU(genome, strlen(genome));
genome = NULL;
// Do transcription and translation
indiv4->do_transcription();
indiv4->do_translation();
}
// ***************************************************************************
// The following commented code allows to print stuff about rnas and proteins
// printf("%"PRId32" rnas and %"PRId32" prots\n", indiv4->get_rna_list()->size(), indiv4->get_protein_list()->size());
// list_node<Rna*>* rna_node = indiv4->get_rna_list()->get_first();
// printf("%"PRId32" rnas and %"PRId32" prots\n", indiv4->rna_list()->size(), indiv4->protein_list()->size());
// list_node<Rna*>* rna_node = indiv4->rna_list()->first();
// while (rna_node != NULL)
// {
// printf("%s rna at pos %"PRId32" (%f, %"PRId32")\n",
// rna_node->get_obj()->get_strand() == LEADING ? "LEADING":"LAGGING",
// rna_node->get_obj()->get_promoter_pos(),
// rna_node->get_obj()->get_basal_level(),
// rna_node->get_obj()->get_transcript_length());
// rna_node->obj()->strand() == LEADING ? "LEADING":"LAGGING",
// rna_node->obj()->promoter_pos(),
// rna_node->obj()->basal_level(),
// rna_node->obj()->transcript_length());
// rna_node = rna_node->get_next();
// rna_node = rna_node->next();
// }
// list_node<Protein*>* protein_node = indiv4->get_protein_list()->get_first();
// list_node<Protein*>* protein_node = indiv4->protein_list()->first();
// while (protein_node != NULL)
// {
// printf("%s protein at pos %"PRId32" (length: %"PRId32", concentr: %f, nb_rnas: %"PRId32")\n",
// protein_node->get_obj()->get_strand() == LEADING ? "LEADING":"LAGGING",
// protein_node->get_obj()->get_shine_dal_pos(),
// protein_node->get_obj()->get_length(),
// protein_node->get_obj()->get_concentration(),
// protein_node->get_obj()->get_rna_list()->size());
// protein_node->obj()->strand() == LEADING ? "LEADING":"LAGGING",
// protein_node->obj()->shine_dal_pos(),
// protein_node->obj()->length(),
// protein_node->obj()->concentration(),
// protein_node->obj()->rna_list()->size());
ExpSetup* exp_s = new ExpSetup(nullptr);
for (int i = 0; i <= 1; i++) {
exp_s->set_fuzzy_flavor(i);
......@@ -258,7 +277,7 @@ void IndividualTest::SetUp(void)
// Build indiv2
// Reverse the whole genome
genome = indiv1->get_genetic_unit(0).get_dna()->get_subsequence(0, 0,
genome = indiv1->genetic_unit(0).dna()->subsequence(0, 0,
LAGGING);
indiv2 = new Individual(nullptr, nullptr, nullptr,
std::make_shared<MutationParams>(params_mut), 1.0,
......@@ -291,7 +310,7 @@ void IndividualTest::SetUp(void)
// Build indiv4
genome = indiv3->get_genetic_unit(0).get_dna()->get_subsequence(0, 0,
genome = indiv3->genetic_unit(0).dna()->subsequence(0, 0,
LAGGING);
indiv4 = new Individual(nullptr, nullptr, nullptr,
std::make_shared<MutationParams>(params_mut), 1.0,
......@@ -307,30 +326,30 @@ void IndividualTest::SetUp(void)
// ***************************************************************************
// The following commented code allows to print stuff about rnas and proteins
// printf("%"PRId32" rnas and %"PRId32" prots\n", indiv4->get_rna_list()->size(), indiv4->get_protein_list()->size());
// list_node<Rna*>* rna_node = indiv4->get_rna_list()->get_first();
// printf("%"PRId32" rnas and %"PRId32" prots\n", indiv4->rna_list()->size(), indiv4->protein_list()->size());
// list_node<Rna*>* rna_node = indiv4->rna_list()->first();
// while (rna_node != NULL)
// {
// printf("%s rna at pos %"PRId32" (%f, %"PRId32")\n",
// rna_node->get_obj()->get_strand() == LEADING ? "LEADING":"LAGGING",
// rna_node->get_obj()->get_promoter_pos(),
// rna_node->get_obj()->get_basal_level(),
// rna_node->get_obj()->get_transcript_length());
// rna_node->obj()->strand() == LEADING ? "LEADING":"LAGGING",
// rna_node->obj()->promoter_pos(),
// rna_node->obj()->basal_level(),
// rna_node->obj()->transcript_length());
// rna_node = rna_node->get_next();
// rna_node = rna_node->next();
// }
// list_node<Protein*>* protein_node = indiv4->get_protein_list()->get_first();
// list_node<Protein*>* protein_node = indiv4->protein_list()->first();
// while (protein_node != NULL)
// {
// printf("%s protein at pos %"PRId32" (length: %"PRId32", concentr: %f, nb_rnas: %"PRId32")\n",
// protein_node->get_obj()->get_strand() == LEADING ? "LEADING":"LAGGING",
// protein_node->get_obj()->get_shine_dal_pos(),
// protein_node->get_obj()->get_length(),
// protein_node->get_obj()->get_concentration(),
// protein_node->get_obj()->get_rna_list()->size());
// protein_node->obj()->strand() == LEADING ? "LEADING":"LAGGING",
// protein_node->obj()->shine_dal_pos(),
// protein_node->obj()->length(),
// protein_node->obj()->concentration(),
// protein_node->obj()->rna_list()->size());
// protein_node = protein_node->get_next();
// protein_node = protein_node->next();
// }
}
}
......@@ -351,115 +370,115 @@ TEST_F(IndividualTest, TestIndiv1)
// "right" means those values we have computed by hand
// Check genome size
EXPECT_EQ(109, indiv1->get_amount_of_dna());
EXPECT_EQ(109, indiv1->get_genetic_unit_seq_length(0));
EXPECT_EQ(109, indiv1->amount_of_dna());
EXPECT_EQ(109, indiv1->genetic_unit_seq_length(0));
// Check RNA list
list<const Rna*> rna_list = indiv1->get_rna_list();
list<const Rna*> rna_list = indiv1->rna_list();
EXPECT_EQ(2, rna_list.size());
const Rna* rna = rna_list.front();
EXPECT_EQ(LEADING, rna->get_strand());
EXPECT_EQ(4, rna->get_promoter_pos());
EXPECT_FLOAT_EQ(0.6, rna->get_basal_level());
EXPECT_EQ(50, rna->get_transcript_length());
EXPECT_EQ(LEADING, rna->strand());
EXPECT_EQ(4, rna->promoter_pos());
EXPECT_FLOAT_EQ(0.6, rna->basal_level());
EXPECT_EQ(50, rna->transcript_length());
rna = rna_list.back();
EXPECT_EQ(LEADING, rna->get_strand());
EXPECT_EQ(81, rna->get_promoter_pos());
EXPECT_FLOAT_EQ(0.8, rna->get_basal_level());
EXPECT_EQ(82, rna->get_transcript_length());
EXPECT_EQ(LEADING, rna->strand());
EXPECT_EQ(81, rna->promoter_pos());
EXPECT_FLOAT_EQ(0.8, rna->basal_level());
EXPECT_EQ(82, rna->transcript_length());
// Check protein list
list<Protein*> prot_list = indiv1->get_protein_list();
list<Protein*> prot_list = indiv1->protein_list();
EXPECT_EQ(1, prot_list.size());
Protein* prot = prot_list.front();
EXPECT_EQ(LEADING, prot->get_strand());
EXPECT_EQ(31, prot->get_shine_dal_pos());
EXPECT_EQ(4, prot->get_length());
EXPECT_FLOAT_EQ(1.4, prot->get_concentration());
EXPECT_EQ(2, prot->get_rna_list().size());
EXPECT_EQ(LEADING, prot->strand());
EXPECT_EQ(31, prot->shine_dal_pos());
EXPECT_EQ(4, prot->length());
EXPECT_FLOAT_EQ(1.4, prot->concentration());
EXPECT_EQ(2, prot->rna_list().size());
}
TEST_F(IndividualTest, TestIndiv2)
{
// Check genome size
EXPECT_EQ(109, indiv2->get_amount_of_dna());
EXPECT_EQ(109, indiv2->get_genetic_unit_seq_length(0));
EXPECT_EQ(109, indiv2->amount_of_dna());
EXPECT_EQ(109, indiv2->genetic_unit_seq_length(0));
// Check RNA list
list<const Rna*> rna_list = indiv2->get_rna_list();
list<const Rna*> rna_list = indiv2->rna_list();
EXPECT_EQ(2, rna_list.size());
const Rna* rna = rna_list.front();
EXPECT_EQ(LAGGING, rna->get_strand());
EXPECT_EQ(104, rna->get_promoter_pos());
EXPECT_FLOAT_EQ(0.6, rna->get_basal_level());
EXPECT_EQ(50, rna->get_transcript_length());
EXPECT_EQ(LAGGING, rna->strand());
EXPECT_EQ(104, rna->promoter_pos());
EXPECT_FLOAT_EQ(0.6, rna->basal_level());
EXPECT_EQ(50, rna->transcript_length());
rna = rna_list.back();
EXPECT_EQ(LAGGING, rna->get_strand());
EXPECT_EQ(27, rna->get_promoter_pos());
EXPECT_FLOAT_EQ(0.8, rna->get_basal_level());
EXPECT_EQ(82, rna->get_transcript_length());
EXPECT_EQ(LAGGING, rna->strand());
EXPECT_EQ(27, rna->promoter_pos());
EXPECT_FLOAT_EQ(0.8, rna->basal_level());
EXPECT_EQ(82, rna->transcript_length());
// Check protein list
list<Protein*> prot_list = indiv2->get_protein_list();
list<Protein*> prot_list = indiv2->protein_list();
EXPECT_EQ(1, prot_list.size());
Protein* prot = prot_list.front();
EXPECT_EQ(LAGGING, prot->get_strand());
EXPECT_EQ(77, prot->get_shine_dal_pos());
EXPECT_EQ(4, prot->get_length());
EXPECT_FLOAT_EQ(1.4, prot->get_concentration());
EXPECT_EQ(2, prot->get_rna_list().size());
EXPECT_EQ(LAGGING, prot->strand());
EXPECT_EQ(77, prot->shine_dal_pos());
EXPECT_EQ(4, prot->length());
EXPECT_FLOAT_EQ(1.4, prot->concentration());
EXPECT_EQ(2, prot->rna_list().size());
}
TEST_F(IndividualTest, TestIndiv3)
{
// Check genome size
EXPECT_EQ(81, indiv3->get_amount_of_dna());
EXPECT_EQ(81, indiv3->get_genetic_unit_seq_length(0));
EXPECT_EQ(81, indiv3->amount_of_dna());
EXPECT_EQ(81, indiv3->genetic_unit_seq_length(0));
// Check RNA list
list<const Rna*> rna_list = indiv3->get_rna_list();
list<const Rna*> rna_list = indiv3->rna_list();
EXPECT_EQ(1, rna_list.size());
const Rna* rna = rna_list.front();
EXPECT_EQ(LEADING, rna->get_strand());
EXPECT_EQ(54, rna->get_promoter_pos());
EXPECT_FLOAT_EQ(0.8, rna->get_basal_level());
EXPECT_EQ(42, rna->get_transcript_length());
EXPECT_EQ(LEADING, rna->strand());
EXPECT_EQ(54, rna->promoter_pos());
EXPECT_FLOAT_EQ(0.8, rna->basal_level());
EXPECT_EQ(42, rna->transcript_length());
// Check protein list
list<Protein*> prot_list = indiv3->get_protein_list();
list<Protein*> prot_list = indiv3->protein_list();
EXPECT_EQ(1, prot_list.size());
Protein* prot = prot_list.front();
EXPECT_EQ(LEADING, prot->get_strand());
EXPECT_EQ(4, prot->get_shine_dal_pos());
EXPECT_EQ(4, prot->get_length());
EXPECT_FLOAT_EQ(0.8, prot->get_concentration());
EXPECT_EQ(1, prot->get_rna_list().size());
EXPECT_EQ(LEADING, prot->strand());
EXPECT_EQ(4, prot->shine_dal_pos());
EXPECT_EQ(4, prot->length());
EXPECT_FLOAT_EQ(0.8, prot->concentration());
EXPECT_EQ(1, prot->rna_list().size());
}
TEST_F(IndividualTest, TestIndiv4)
{
// Check genome size
EXPECT_EQ(81, indiv4->get_amount_of_dna());
EXPECT_EQ(81, indiv4->get_genetic_unit_seq_length(0));
EXPECT_EQ(81, indiv4->amount_of_dna());
EXPECT_EQ(81, indiv4->genetic_unit_seq_length(0));
// Check RNA list
list<const Rna*> rna_list = indiv4->get_rna_list();
list<const Rna*> rna_list = indiv4->rna_list();
EXPECT_EQ(1, rna_list.size());
const Rna* rna = rna_list.front();
EXPECT_EQ(LAGGING, rna->get_strand());
EXPECT_EQ(26, rna->get_promoter_pos());
EXPECT_FLOAT_EQ(0.8, rna->get_basal_level());
EXPECT_EQ(42, rna->get_transcript_length());
EXPECT_EQ(LAGGING, rna->strand());
EXPECT_EQ(26, rna->promoter_pos());
EXPECT_FLOAT_EQ(0.8, rna->basal_level());
EXPECT_EQ(42, rna->transcript_length());
// Check protein list
list<Protein*> prot_list = indiv4->get_protein_list();
list<Protein*> prot_list = indiv4->protein_list();
EXPECT_EQ(1, prot_list.size());
Protein* prot = prot_list.front();
EXPECT_EQ(LAGGING, prot->get_strand());
EXPECT_EQ(76, prot->get_shine_dal_pos());
EXPECT_EQ(4, prot->get_length());
EXPECT_FLOAT_EQ(0.8, prot->get_concentration());
EXPECT_EQ(1, prot->get_rna_list().size());
EXPECT_EQ(LAGGING, prot->strand());
EXPECT_EQ(76, prot->shine_dal_pos());
EXPECT_EQ(4, prot->length());
EXPECT_FLOAT_EQ(0.8, prot->concentration());
EXPECT_EQ(1, prot->rna_list().size());
}
// ===========================================================================
......
src/tests/TranscriptionTranslationTest.cpp
View file @ 7ffdc5fd
......@@ -114,29 +114,29 @@ list<protein_line> analyse_gu(const GeneticUnit& gen_unit, int32_t gen_unit_numb
{
list<protein_line> proteins;
Promoters2Strands llrnas = gen_unit.get_rna_list();
Promoters2Strands llrnas = gen_unit.rna_list();
for(auto lrnas : llrnas) {
for (auto rna : lrnas) {
for (auto prot : rna.get_transcribed_proteins()) {
double mean = prot->get_mean();
for (auto prot : rna.transcribed_proteins()) {
double mean = prot->mean();
int nfeat = -1;
protein_line prot_line;
prot_line.id = gen_unit.get_indiv()->get_id();
prot_line.id = gen_unit.indiv()->id();
prot_line.c_or_p = gen_unit_number != 0 ? "PLASMID" : "CHROM";
prot_line.strand = prot->get_strand() == LEADING ? "LEADING" : "LAGGING";
prot_line.pos = prot->get_first_translated_pos();
prot_line.len = prot->get_length();
prot_line.lpos = prot->get_last_translated_pos();
prot_line.sequence = prot->get_AA_sequence('_');
prot_line.strand = prot->strand() == LEADING ? "LEADING" : "LAGGING";
prot_line.pos = prot->first_translated_pos();
prot_line.len = prot->length();
prot_line.lpos = prot->last_translated_pos();
prot_line.sequence = prot->AA_sequence('_');
prot_line.m = mean;
prot_line.w = prot->get_width();
prot_line.h = prot->get_height();
prot_line.c = prot->get_concentration();
prot_line.w = prot->width();
prot_line.h = prot->height();
prot_line.c = prot->concentration();
prot_line.f = nfeat;
prot_line.prom_pos = rna.get_promoter_pos();
prot_line.rna_len = rna.get_transcript_length();
prot_line.basal_level = rna.get_basal_level();
prot_line.prom_pos = rna.promoter_pos();
prot_line.rna_len = rna.transcript_length();
prot_line.basal_level = rna.basal_level();
proteins.push_back(prot_line);
}
......@@ -170,8 +170,8 @@ void expect_equal(const list<protein_line> expected_proteins,
}
}
void TranscriptionTranslationTest::check_genome(const string& dir, int generation) {
void TranscriptionTranslationTest::check_genome(const string& dir,
int generation) {
ExpSetup* exp_s = new ExpSetup(nullptr);
for (int i = 0; i <= 1; i++) {
exp_s->set_fuzzy_flavor(i);
......@@ -209,7 +209,7 @@ void TranscriptionTranslationTest::check_genome(const string& dir, int generatio
indiv->do_transcription();
indiv->do_translation();
auto actual_proteins = analyse_gu(indiv->get_genetic_unit(0), 0);
auto actual_proteins = analyse_gu(indiv->genetic_unit(0), 0);
// Read the expected results
......@@ -231,7 +231,6 @@ void TranscriptionTranslationTest::check_genome(const string& dir, int generatio
// (very dirty) Get rid of last line (added twice)
expected_proteins.pop_back();
// cout << "*************** EXPECTED ********************" << endl;
// for (auto prot_line : expected_proteins) {
// cout << prot_line << endl;
......@@ -243,20 +242,6 @@ void TranscriptionTranslationTest::check_genome(const string& dir, int generatio
expect_equal(expected_proteins, actual_proteins);
}
// (very dirty) Get rid of last line (added twice)
expected_proteins.pop_back();
// cout << "*************** EXPECTED ********************" << endl;
// for (auto prot_line : expected_proteins) {
// cout << prot_line << endl;
// }
// cout << "**************** ACTUAL *********************" << endl;
// for (auto prot_line : actual_proteins) {
// cout << prot_line << endl;
// }
expect_equal(expected_proteins, actual_proteins);
}
TEST_F(TranscriptionTranslationTest, TestIndivVirus6) {
......