Adding support for undetermined AA from codon containing N

cpp/Vihgo.cpp

@@ -550,7 +550,7 @@ void Vihgo::writeCodonFile( string incType )
 				{
 					string ref = *myIterD;
 					string refCodon = this->referenceCodon[myIterA->first][myIterB->first];
-					string refAA = this->codonTable[refCodon];
+					string refAA = this->getAAfromCodon(refCodon);
 					int position = this->refProtCorresp[ref][myIterB->first];
 					double depth = 0.0;
 
@@ -563,7 +563,7 @@ void Vihgo::writeCodonFile( string incType )
 					for( myIterC = (myIterB->second).begin() ; myIterC != (myIterB->second).end() ; myIterC ++ )
 					{
 						string altCodon = myIterC->first;
-						string altAA = this->codonTable[altCodon];
+						string altAA = this->getAAfromCodon(altCodon);
 						double altCount = myIterC->second;
 						double allelicRatio = ( altCount / depth ) * 100;
 
@@ -610,6 +610,7 @@ void Vihgo::writeCodonFile( string incType )
 }
 
 
+
 // Method to get annotation for variation
 // Actually has been deactivated yet
 // NEED A REWAMP FOR ADAPTATION TO ALL VIRUS; VCF FORMAT
@@ -1058,14 +1059,14 @@ int Vihgo::variantCalling()
 				double altNumber = 0.0;
 
 				string refCodon = this->referenceCodon[myIterA->first][myIterB->first];
-				string aa = codonTable[refCodon];
+				string aa = this->getAAfromCodon(refCodon);
 
 
 				// cerr << "Analysing codon " <<  myIterA->first << ":" << myIterB->first << " ; Ref is : " << refCodon  << " : " << aa << endl;
 				for( myIterC = (myIterB->second).begin() ; myIterC != (myIterB->second).end() ; myIterC ++ )
 				{
 					string currentCodon = myIterC->first;
-					string currentAa = codonTable[currentCodon];
+					string currentAa = this->getAAfromCodon(currentCodon);
 
 					// cerr << "\tin control : " << currentCodon << " : " << currentAa << " : " << myIterC->second << endl;
 
@@ -1107,14 +1108,14 @@ int Vihgo::variantCalling()
 			double altNumber = 0.0;
 
 			string refCodon = this->referenceCodon[myIterA->first][myIterB->first];
-			string aa = codonTable[refCodon];
+			string aa = this->getAAfromCodon(refCodon);
 
 
 			// cerr << "Analysing codon " <<  myIterA->first << ":" << myIterB->first << " ; Ref is : " << refCodon  << " : " << aa << endl;
 			for( myIterC = (myIterB->second).begin() ; myIterC != (myIterB->second).end() ; myIterC ++ )
 			{
 				string currentCodon = myIterC->first;
-				string currentAa = codonTable[currentCodon];
+				string currentAa = this->getAAfromCodon(currentCodon);
 
 				// cerr << "\tin sample : " << currentCodon << " : " << currentAa << " : " << myIterC->second << endl;
 
@@ -1240,7 +1241,7 @@ int Vihgo::variantCalling()
 	// result file
 	ofstream outresultStream;
 	outresultStream.open( (this->outputFile).c_str() , ios::out);
-	outresultStream <<  "#gene\tchromosome\tgPosition\tAAPosition\treferenceCodon\treferenceAA\tvariantCodon\tvariantAA\ttotalDepth\tAlternativeAlleleratio\tFETp-value\tFETbgp-value\tAnnotation" << endl;
+	outresultStream <<  "#chromosome\tgPosition\tgene\tAAPosition\treferenceCodon\treferenceAA\tvariantCodon\tvariantAA\ttotalDepth\tAlternativeAlleleratio\tFETp-value\tFETbgp-value\tAnnotation" << endl;
 
 	// test each position
 	cerr << "Variant calling in progress ... " << endl;
@@ -1282,11 +1283,11 @@ int Vihgo::variantCalling()
 
 				// cerr << "\t\tgetting ref informations ..." << endl;
 				string refCodon = this->referenceCodon[myIterD->first][myIterE->first];
-				string aa = codonTable[refCodon];
+				string aa = this->getAAfromCodon(refCodon);
 
 				// cerr << "\t\tgetting sample informations ..." << endl;
 				string currentCodon = myIterF->first;
-				string currentAa = codonTable[currentCodon];
+				string currentAa = this->getAAfromCodon(currentCodon);
 
 				// disengage variant calling if ref aa = alt aa (synonymous)
 				if (aa == currentAa)
@@ -1923,7 +1924,18 @@ char Vihgo::getReferenceBase( string incChr , int incPos )
 	}
 }
 
-
+// method to get AA from a codon
+string Vihgo::getAAfromCodon( string incCodon )
+{
+	if ( this->codonTable.count( incCodon ) > 0)
+	{
+		return this->codonTable[incCodon];
+	}
+	else
+	{
+		return "?";
+	}
+}
 
 
 vector<string> Vihgo::getAnnotationFiles()

cpp/Vihgo.h

@@ -73,6 +73,8 @@ class Vihgo
 		char getBase( int );
 		char getReferenceBase( std::string , int );
 
+
+
 		// map of codon
 		std::map<std::string , std::string > codonTable;