Added initial Avida code (computes 0-dimensional persistence of the hamming distance function)
authorDmitriy Morozov <morozov@cs.duke.edu>
Thu, 27 Dec 2007 00:00:00 -0500
changeset 43 0f970f454094
parent 40 122e4a1fa117
child 44 0bf12743641f
Added initial Avida code (computes 0-dimensional persistence of the hamming distance function)
examples/CMakeLists.txt
examples/cech-complex/cech-complex.cpp
examples/fitness/CMakeLists.txt
examples/fitness/avida-distance.cpp
examples/fitness/avida-population-detail.h
examples/fitness/avida-population-detail.hpp
include/topology/simplex.h
--- a/examples/CMakeLists.txt	Sun Nov 25 09:12:36 2007 -0500
+++ b/examples/CMakeLists.txt	Thu Dec 27 00:00:00 2007 -0500
@@ -1,5 +1,6 @@
 add_subdirectory			(alphashapes)
 add_subdirectory			(ar-vineyard)
 add_subdirectory			(cech-complex)
+add_subdirectory			(fitness)
 add_subdirectory			(grid)
 add_subdirectory			(triangle)
--- a/examples/cech-complex/cech-complex.cpp	Sun Nov 25 09:12:36 2007 -0500
+++ b/examples/cech-complex/cech-complex.cpp	Thu Dec 27 00:00:00 2007 -0500
@@ -14,18 +14,7 @@
 typedef			SimplexWithValue<PointIndex>							Simplex;
 typedef 		std::vector<Simplex> 									SimplexVector;
 typedef 		Filtration<Simplex>										CechFiltration;
-
-class DimensionValueComparison
-{
-	public:
-		bool	operator()(const Simplex& s1, const Simplex& s2) const
-		{
-			if (s1.dimension() == s2.dimension())
-				return s1.get_value() < s2.get_value();
-			else
-				return s1.dimension() < s2.dimension();
-		}
-};
+typedef			DimensionValueComparison<Simplex>						DimValComparison;
 
 int choose(int n, int k)
 {
@@ -109,7 +98,7 @@
 			add_simplices(sv, i, points);
 		std::cout << "Size of SimplexVector: " << sv.size() << std::endl;
 			
-		std::sort(sv.begin(), sv.end(), DimensionValueComparison());
+		std::sort(sv.begin(), sv.end(), DimValComparison());
 		
 		for (SimplexVector::const_iterator cur = sv.begin(); cur != sv.end(); ++cur)
 			cf.append(*cur);
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/examples/fitness/CMakeLists.txt	Thu Dec 27 00:00:00 2007 -0500
@@ -0,0 +1,7 @@
+set							(targets						
+							 avida)
+							 
+foreach 					(t ${targets})
+	add_executable			(${t} ${t}.cpp)
+	target_link_libraries	(${t} ${libraries})
+endforeach 					(t ${targets})
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/examples/fitness/avida-distance.cpp	Thu Dec 27 00:00:00 2007 -0500
@@ -0,0 +1,94 @@
+#include <iostream>
+#include <vector>
+#include <algorithm>
+#include "avida-population-detail.h"
+
+#include <topology/filtration.h>
+#include <topology/simplex.h>
+
+
+typedef			SimplexWithValue<AvidaOrganismDetail::IDType>		Simplex;
+typedef			std::vector<Simplex>								SimplexVector;
+typedef			Filtration<Simplex>									SimplexFiltration;
+
+int main(int argc, char** argv)
+{
+#ifdef LOGGING
+	rlog::RLogInit(argc, argv);
+	//stdoutLog.subscribeTo(RLOG_CHANNEL("info"));
+#endif
+
+	if (argc < 2)
+	{
+		std::cout << "USAGE: avida FILENAME" << std::endl;
+		return 0;
+	}
+
+	AvidaPopulationDetail population(argv[1]);
+	const AvidaPopulationDetail::OrganismVector& organisms = population.get_organisms();
+
+	rInfo("Number of organisms: %d", organisms.size());
+	for (int i = 0; i < population.get_organisms().size(); ++i)
+		rInfo("%d (%s) %f %d %d", organisms[i].id(),
+								  organisms[i].genome().c_str(),
+								  organisms[i].fitness(),
+								  organisms[i].length(),
+								  organisms[i].genome().size());
+
+	// Distance function filtration
+	SimplexVector simplices;
+
+	// Insert edges
+	AvidaOrganismDetail::DistanceType avg_distance = 0;
+	for (AvidaOrganismDetail::CountType i = 0; i < organisms.size(); ++i)
+	{
+		simplices.push_back(0);
+		simplices.back().add(organisms[i].id());
+
+		for (AvidaOrganismDetail::CountType j = i+1; j < organisms.size(); ++j)
+		{
+			avg_distance += organisms[i].genome_distance(organisms[j]);
+			simplices.push_back(Simplex(organisms[i].genome_distance(organisms[j])));
+			simplices.back().add(organisms[i].id());
+			simplices.back().add(organisms[j].id());
+		}
+	}
+	std::sort(simplices.begin(), simplices.end(), DimensionValueComparison<Simplex>());
+	rInfo("Average distance: %f", float(avg_distance)/
+								  ((organisms.size()*organisms.size() - organisms.size())/2));
+
+	SimplexFiltration filtration;
+	for (SimplexVector::const_iterator cur = simplices.begin(); cur != simplices.end(); ++cur)
+	{
+		rInfo("Simplex: %s", tostring(*cur).c_str());
+		filtration.append(*cur);
+	}
+
+	filtration.fill_simplex_index_map();
+	filtration.pair_simplices(false);			// pair simplices without storing trails
+
+	std::cout << "Outputting histogram of death values" << std::endl;
+	typedef std::vector<RealType> DeathVector;
+	DeathVector deaths;
+	for (SimplexFiltration::Index i = filtration.begin(); i != filtration.end(); ++i)
+	{
+		if (i->is_paired())
+			if (i->sign())
+			{
+				AssertMsg(i->dimension() == 0, "Expecting only 0-dimensional diagram");
+				AssertMsg(i->get_value() == 0, "Expecting only 0 birth values in 0-D diagram ");
+				deaths.push_back(i->pair()->get_value());
+			}
+	}
+
+	// Produce histogram
+	std::sort(deaths.begin(), deaths.end());
+	for (DeathVector::iterator cur = deaths.begin(); cur != deaths.end(); )
+	{
+		DeathVector::iterator nw = std::find_if(cur, deaths.end(), 
+												std::bind2nd(std::greater<RealType>(), *cur));
+		std::cout << *cur << "\t" << (nw - cur) << std::endl;
+		cur = nw;
+	}
+	std::cout << "Total: " << deaths.size() + 1;		// +1 for the unpaired
+}
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/examples/fitness/avida-population-detail.h	Thu Dec 27 00:00:00 2007 -0500
@@ -0,0 +1,62 @@
+/**
+ * Author: Dmitriy Morozov
+ * Department of Computer Science, Duke University, 2007
+ */
+
+#ifndef __AVIDA_POPULATION_DETAIL_H__
+#define __AVIDA_POPULATION_DETAIL_H__
+
+#include <string>
+#include <vector>
+
+/**
+ * Stores organism details stored in a single line of a population detail file in data/ 
+ * directory of Avida's output.
+ */
+class AvidaOrganismDetail
+{
+	public:
+		typedef			int							IDType;
+		typedef			unsigned int				DistanceType;
+		typedef			unsigned int				CountType;
+
+						AvidaOrganismDetail(std::string line);
+
+		DistanceType 	genome_distance(const AvidaOrganismDetail& other) const;
+		
+		IDType			id() const									{ return id_; }
+		float			fitness() const								{ return fitness_; }
+		CountType		length() const								{ return genome_length_; }
+		std::string		genome() const								{ return genome_; }
+
+	private:
+		IDType			id_, parent_id_;
+		int				parent_distance_;
+		CountType		num_organisms_alive_, num_organisms_ever_;
+		CountType		genome_length_;
+		float 			merit_, gestation_time_;
+		float 			fitness_;
+		int 			update_born_, update_deactivated_, depth_phylogenetic_tree_;
+		std::string 	genome_;
+};
+
+/**
+ * Stores entire population details (all organisms in a given time step), i.e., stores
+ * an entire population detail file in data/ directory of Avida's ouptut.
+ */
+class AvidaPopulationDetail
+{
+	public:
+		typedef			std::vector<AvidaOrganismDetail>				OrganismVector;
+						AvidaPopulationDetail(std::string filename);
+
+		const OrganismVector& get_organisms() const						{ return organisms_; }
+
+	private:
+		OrganismVector	organisms_;
+};
+
+
+#include "avida-population-detail.hpp"
+
+#endif //__AVIDA_POPULATION_DETAIL_H__
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/examples/fitness/avida-population-detail.hpp	Thu Dec 27 00:00:00 2007 -0500
@@ -0,0 +1,57 @@
+#include <fstream>
+#include <sstream>
+#include "utilities/log.h"
+
+/** AvidaOrganismDetail */
+
+AvidaOrganismDetail::
+AvidaOrganismDetail(std::string line)
+{
+	std::stringstream linestream(line);
+	linestream >> id_ >> parent_id_ >> parent_distance_;
+	linestream >> num_organisms_alive_ >> num_organisms_ever_;
+	linestream >> genome_length_;
+	linestream >> merit_ >> gestation_time_;
+	linestream >> fitness_;
+	linestream >> update_born_ >> update_deactivated_ >> depth_phylogenetic_tree_;
+	linestream >> genome_;
+
+	AssertMsg(genome_length_ == genome_.size(), "Genome must be of given length");
+}
+
+AvidaOrganismDetail::DistanceType
+AvidaOrganismDetail::
+genome_distance(const AvidaOrganismDetail& other) const
+{
+	AssertMsg(genome_.size() == other.genome_.size(), "Currently genome sizes must be the same for distance computation");
+	AssertMsg(genome_length_ == genome_.size(), "Genome length value must match the length of the genome string");
+	AssertMsg(other.genome_length_ == other.genome_.size(), "Genome length value must match the length of the genome string");
+
+	CountType count = 0;
+	for (CountType i = 0; i < genome_.size(); ++i)
+		if (genome_[i] != other.genome_[i])
+			++count;
+
+	return count;
+}
+
+
+/** AvidaPopulationDetail */
+AvidaPopulationDetail::
+AvidaPopulationDetail(std::string filename)
+{
+	std::ifstream infile(filename.c_str());
+	while(infile)
+	{
+		std::string line;
+		std::getline(infile, line);
+		
+		// Skip comments and empty lines
+		char c = '#';
+		std::istringstream linestream(line);
+		linestream >> c;
+		if (c == '#') continue;
+
+		organisms_.push_back(AvidaOrganismDetail(line));
+	}
+}
--- a/include/topology/simplex.h	Sun Nov 25 09:12:36 2007 -0500
+++ b/include/topology/simplex.h	Thu Dec 27 00:00:00 2007 -0500
@@ -165,6 +165,21 @@
 		void 					serialize(Archive& ar, version_type );
 };
 
+template<class Simplex_>
+class DimensionValueComparison
+{
+	public:
+		typedef					Simplex_									Simplex;
+
+		bool					operator()(const Simplex& s1, const Simplex& s2) const
+		{
+			if (s1.dimension() == s2.dimension())
+				return s1.get_value() < s2.get_value();
+			else
+				return s1.dimension() < s2.dimension();
+		}
+};
+
 #include "simplex.hpp"
 
 #endif // __SIMPLEX_H__