[Avida-SVN] r1268 - in trunk: Avida.xcodeproj source/analyze
matt at myxo.css.msu.edu
matt at myxo.css.msu.edu
Mon Feb 12 21:51:12 PST 2007
Author: matt
Date: 2007-02-13 00:51:12 -0500 (Tue, 13 Feb 2007)
New Revision: 1268
Modified:
trunk/Avida.xcodeproj/project.pbxproj
trunk/source/analyze/cAnalyze.cc
trunk/source/analyze/cAnalyze.h
Log:
Removed previous, faulty versions of epistasis analysis.
Modified new AnalyzeEpistasis function to keep single
site probabilities constant during each iterative update
cycle.
Modified: trunk/Avida.xcodeproj/project.pbxproj
===================================================================
--- trunk/Avida.xcodeproj/project.pbxproj 2007-02-10 21:37:24 UTC (rev 1267)
+++ trunk/Avida.xcodeproj/project.pbxproj 2007-02-13 05:51:12 UTC (rev 1268)
@@ -810,7 +810,7 @@
DCC315CF076253A5008F7A48 /* Makefile */ = {isa = PBXFileReference; fileEncoding = 30; lastKnownFileType = sourcecode.make; path = Makefile; sourceTree = "<group>"; };
DCC315D0076253A5008F7A48 /* task_event_gen.cc */ = {isa = PBXFileReference; fileEncoding = 30; lastKnownFileType = sourcecode.cpp.cpp; path = task_event_gen.cc; sourceTree = "<group>"; };
DCC315D1076253A5008F7A48 /* task_event_gen.old.cc */ = {isa = PBXFileReference; fileEncoding = 30; lastKnownFileType = sourcecode.cpp.cpp; path = task_event_gen.old.cc; sourceTree = "<group>"; };
- DCC3164D07626CF3008F7A48 /* avida */ = {isa = PBXFileReference; includeInIndex = 0; lastKnownFileType = "compiled.mach-o.executable"; path = avida; sourceTree = BUILT_PRODUCTS_DIR; };
+ DCC3164D07626CF3008F7A48 /* avida */ = {isa = PBXFileReference; explicitFileType = "compiled.mach-o.executable"; includeInIndex = 0; path = avida; sourceTree = BUILT_PRODUCTS_DIR; };
/* End PBXFileReference section */
/* Begin PBXFrameworksBuildPhase section */
Modified: trunk/source/analyze/cAnalyze.cc
===================================================================
--- trunk/source/analyze/cAnalyze.cc 2007-02-10 21:37:24 UTC (rev 1267)
+++ trunk/source/analyze/cAnalyze.cc 2007-02-13 05:51:12 UTC (rev 1268)
@@ -569,128 +569,9 @@
return entropy;
}
-//@ MRR @@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
-tMatrix< double > cAnalyze::AnalyzeEntropyPairs(cAnalyzeGenotype * genotype, double mu)
-{
-
- double entropy = 0.0;
- genotype->Recalculate(m_ctx, m_testcpu);
- // Calculate the stats for the genotype we're working with ...
- const int num_insts = inst_set.GetSize();
- const int num_lines = genotype->GetLength();
- const cGenome & base_genome = genotype->GetGenome();
- cGenome mod_genome(base_genome);
- double base_fitness = genotype->GetFitness();
- cout << num_lines << endl;
- tMatrix< double > pairwiseEntropy(num_lines, num_lines);
- for (int i=0; i<num_lines; i++)
- for (int j=-0; j<num_lines; j++)
- pairwiseEntropy[i][j] = 0.0;
-
- cout << pairwiseEntropy.GetNumRows() << endl;
-
- // If the fitness is 0, return empty matrix
-
- if (genotype->GetFitness() == 0) {
- return pairwiseEntropy;
- }
-
-
- tMatrix< double > test_fitness(num_insts,num_insts);
- tMatrix< double > prob(num_insts,num_insts);
-
- //Pairwise mutations; the diagonal of the matrix will be the information
- //stored by that site alone
- for (int line_1 = 0; line_1 < num_lines; line_1++){
- for (int line_2 = line_1; line_2 < num_lines; line_2++) {
-
- cerr << "[ " << line_1 << ", " << line_2 << " ]" << endl;
-
- int cur_inst_1 = base_genome[line_1].GetOp();
- int cur_inst_2 = base_genome[line_2].GetOp();
-
- // Test fitness of each mutant.
- for (int mod_inst_1 = 0; mod_inst_1 < num_insts; mod_inst_1++){
- for (int mod_inst_2 = 0; mod_inst_2 < num_insts; mod_inst_2++) {
- mod_genome[line_1].SetOp(mod_inst_1);
- mod_genome[line_2].SetOp(mod_inst_2);
- cAnalyzeGenotype test_genotype(m_world, mod_genome, inst_set);
- test_genotype.Recalculate(m_ctx, m_testcpu);
- // Adjust fitness ...
- if (test_genotype.GetFitness() <= base_fitness) {
- test_fitness[mod_inst_1][mod_inst_2] = test_genotype.GetFitness();
- } else {
- test_fitness[mod_inst_1][mod_inst_2] = base_fitness;
- }
- }
- }
-
- // Calculate probabilities at mut-sel balance
- double w_bar = 1;
-
- // Normalize fitness values
- double maxFitness = 0.0;
- for(int i=0; i<num_insts; i++) {
- for (int j = 0; j < num_insts; j++){
- if(test_fitness[i][j] > maxFitness) {
- maxFitness = test_fitness[i][j];
- }
- }
- }
-
-
- for(int i=0; i<num_insts; i++) {
- for (int j=0; j<num_insts; j++){
- test_fitness[i][j] /= maxFitness;
- }
- }
-
-
- while(1) {
- double sum = 0.0;
- for (int mod_inst_1 = 0; mod_inst_1 < num_insts; mod_inst_1++) {
- for (int mod_inst_2 = 0; mod_inst_2 < num_insts; mod_inst_2++){
-
- prob[mod_inst_1][mod_inst_2] =
- (mu * w_bar) /
- ((double)num_insts *
- (w_bar + test_fitness[mod_inst_1][mod_inst_2]
- * mu - test_fitness[mod_inst_1][mod_inst_2]));
- sum = sum + prob[mod_inst_1][mod_inst_2];
- }
- }
- if ((sum-1.0)*(sum-1.0) <= 0.0001)
- break;
- else
- w_bar = w_bar - 0.000001;
- }
-
- // Calculate entropy ...
- double this_entropy = 0.0;
- for (int i = 0; i < num_insts; i++){
- for (int j = 0; j < num_insts; j++) {
- this_entropy += prob[i][j] *
- log((double) 1.0/prob[i][j]) / log ((double) num_insts);
- }
- }
- entropy += this_entropy;
- pairwiseEntropy[line_1][line_2] = this_entropy;
-
- // Reset the mod_genome back to the original sequence.
- mod_genome[line_1].SetOp(cur_inst_1);
- mod_genome[line_2].SetOp(cur_inst_2);
-
- }
- } //End Loops
- return pairwiseEntropy;
-}
-
-
-
-
double cAnalyze::AnalyzeEntropyGivenParent(cAnalyzeGenotype * genotype,
cAnalyzeGenotype * parent, double mut_rate)
{
@@ -3429,64 +3310,7 @@
-//@ MRR @@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
-void cAnalyze::CommandPairwiseEntropy(cString cur_string)
-{
- if (m_world->GetVerbosity() >= VERBOSE_ON)
- cout << "Finding pairwise entropy on batch " << cur_batch << endl;
- else
- cout << "Finding pairwise entropy..." << endl;
-
- cout << "@MRR-> This command is being tested!" << endl;
- cString directory = PopDirectory(cur_string, "pairwise_data/");
- if (m_world->GetVerbosity() >= VERBOSE_ON)
- cout << "\tUsing directory: " << directory << endl;
- double mu = cur_string.PopWord().AsDouble();
- if (m_world->GetVerbosity() >= VERBOSE_ON)
- cout << "\tUsing mu=" << mu << endl;
-
- tListIterator<cAnalyzeGenotype> batch_it(batch[cur_batch].List());
- cAnalyzeGenotype * genotype = batch_it.Next();
-
- cout << genotype->GetName() << endl;
-
- while(genotype != NULL)
- {
- cString genName = genotype->GetName();
-
- if (m_world->GetVerbosity() >= VERBOSE_ON)
- cout << "\t...on genotype " << genName << endl;
-
- cString filename;
- filename.Set("%spairdata.%s.dat", static_cast<const char*>(directory),
- static_cast<const char*>(genName));
-
- // @DMB -- ofstream& fp = m_world->GetDataFileOFStream(filename);
-
- if (m_world->GetVerbosity() >= VERBOSE_ON)
- cout << "\t\t...with filename: " << filename << endl;
-
- cout << "# Pairwise Entropy Information" << endl;
-
- tMatrix<double> pairdata = AnalyzeEntropyPairs(genotype, mu);
-
- cout << pairdata.GetNumRows() << endl;
-
- for (int i=0; i < pairdata.GetNumRows(); i++){
- for (int j=0; j < pairdata.GetNumCols(); j++)
- cout << pairdata[i][j] << " ";
- cout << endl;
- }
- m_world->GetDataFileManager().Remove(filename);
- genotype = batch_it.Next();
- }
-}
-
-
-
-
-
// This command will take the current batch and analyze how well organisms
// cross-over with each other, both across the population and between mates.
@@ -6794,17 +6618,21 @@
while(1) {
double sum = 0.0;
+ double Pa = 0.0;
+ double Pb = 0.0;
+ //Collect our Pa & Pb data before we modify it
+ for (int j = 0; j < num_insts; j++){
+ for (int k = 0; k < num_insts; k++){
+ Pa += prob[k][j] * test_fitness[k][j];
+ Pb += prob[j][k] * test_fitness[j][k];
+ }
+ }
+
for (int mod_A = 0; mod_A < num_insts; mod_A++) {
for (int mod_B = 0; mod_B < num_insts; mod_B++){
- double Pa = 0.0;
- double Pb = 0.0;
- for (int k = 0; k < num_insts; k++)
- {
- Pa += prob[k][mod_B] * test_fitness[k][mod_B];
- Pb += prob[mod_A][k] * test_fitness[mod_A][k];
- }
+
double Pab = ( (prob[mod_A][mod_B] * test_fitness[mod_A][mod_B] * mut_elitist)
+ Pa * mut_single_in / num_insts + Pb * mut_single_in / num_insts ) / w_bar
+ mut_double_in / (num_insts*num_insts);
@@ -7857,7 +7685,6 @@
AddLibraryDef("ANALYZE_KNOCKOUTS", &cAnalyze::AnalyzeKnockouts);
AddLibraryDef("ANALYZE_POP_COMPLEXITY", &cAnalyze::AnalyzePopComplexity);
AddLibraryDef("MAP_DEPTH", &cAnalyze::CommandMapDepth);
- AddLibraryDef("PAIRWISE_ENTROPY", &cAnalyze::CommandPairwiseEntropy); //@MRR
// Population comparison commands...
AddLibraryDef("HAMMING", &cAnalyze::CommandHamming);
Modified: trunk/source/analyze/cAnalyze.h
===================================================================
--- trunk/source/analyze/cAnalyze.h 2007-02-10 21:37:24 UTC (rev 1267)
+++ trunk/source/analyze/cAnalyze.h 2007-02-13 05:51:12 UTC (rev 1268)
@@ -171,7 +171,6 @@
void CommandAnalyzeModularity(cString cur_string);
void CommandMapMutations(cString cur_string);
void CommandMapDepth(cString cur_string);
- void CommandPairwiseEntropy(cString cur_string);
// Population Comparison Commands...
void CommandHamming(cString cur_string);
@@ -240,9 +239,6 @@
cAnalyzeGenotype * genotype2,
double mut_rate);
- //Calculate covarying information between pairs of sites
- tMatrix<double> AnalyzeEntropyPairs(cAnalyzeGenotype * genotype,
- double mut_rate);
// Flow Control...
More information about the Avida-cvs
mailing list