
% 
% \subsubsection{Processing annotated genomes from NCBI} 
% 
% The objective here is to generate sets  of genes from each  genome so that
% genes are organized  without any duplication.  The input  is a list of
% chloroplasts annotated genomes downloaded from the NCBI website. All genomes,
% which consist in collections of
% protein  coding sequences~\cite{parra2007cegma,RDOGMA},
% are stored as \textit{fasta} files.   To be able to build the set of
% core    genes,     we    need    to     preprocess    these    genomes
% using  \textit{BioPython}  package \cite{chapman2000biopython}.   
% % 
% First of all, we starts by  converting  each  genome from  \textit{fasta}  file  to
% GenVision \cite{geneVision}  format from DNASTAR. Each  genome is thus
% converted in a list of genes, containing both gene names and occurrences. 
% Gene
% name duplication can be accumulated during the treatment of a genome.
% % 
% These  duplication  come   from  gene  fragments  (\emph{e.g.}   gene
% fragments treated  with NCBI) and  from chloroplast DNA  sequences. To
% ensure that  all the  duplication are removed,  each list of  gene is
% translated  into a  set of  genes.  
% % Note that  NCBI genome  annotation
% produces genes except \textit{Ribosomal (rRNA)} genes.
% 

%\subsubsection{Processing genomes annotated by DOGMA}

% Protein  coding genes are  identified in an  input genome
%using sequence similarity of genes  in DOGMA database.  In addition in
%comparison   with   NCBI    annotation   tool,   
% It   can   detect
% both \textit{transfer RNAs (tRNA)} and \textit{ribosomal RNAs (rRNA)}. 
%  
% The DOGMA  annotation process  is divided into  two tasks.   First, we
% manually annotate chloroplast genomes using DOGMA web tool. The output
% of this step is supposed to  be a collection of coding genes files for
% each genome, organized in GeneVision file. The second task is to solve
% the  gene   duplication  problem  and   therefore  we  have   used  two
% methods. 

\begin{figure}[H]
\centering
\includegraphics[scale=0.3]{gensim}
\caption{Part of the implementation of the second method, compaire the common genes from NCBI and DOGMA.}
\label{Meth2:gensim}
\end{figure}

\color{red}The second approach in this paper is an enhasement of the ICM\textit{Intersection Core Matrix} proposed in \cite{Alkindy2014} by considering gene names to find core genome. Based on gene names spelling, When they realizing simple homogenization of names provided by NCBI, they miss core genes which have slightly different name formats. \color{black} To enlarge the size of the core genome, to be as close as possible to the true natural one, we propose to integrate a similarity distance on gene names. Each similarity will be computed between a name from DOGMA, which operates as a reference here, and a name from NCBI as shown in figure~\ref{Meth2:gensim}. 

The proposed distance  is the Levenshtein one, which is close to the Needleman-Wunsch, except that gap opening and extension penalties are equal. The same name is then set to sequences whose NCBI names are close according to this edit distance.

The risk, by doing so, is to merge genes that are different but whose names are similar (for instance, ND4 and ND4L are two different mitochondrial genes but with similar names). The solution is thus to compare, in a second stage, the similarity of DNA sequences too (with a Needleman-Wunsch global alignment), and to simply ignore the gene if this similarity is below a given threshold. 

By doing so, the second approach is designed, which takes the fundamental idea contained in the annotation-based approaches in the previous work. Remark that this approach is simply a deeper processing of the naming stage in the second approach in \cite{Alkindy2014}, the other stages being identical.

The DNA similarity computation raises another problem in the case of DOGMA:
contrary to what happens with gene features in NCBI, genes predicted by
DOGMA may be fragmented in several parts. Such genes are signaled
in the GeneVision file produced by DOGMA, as each fragment is in this file
and with the same gene name. A gene whose name is present at least twice
in the file is thus either a duplicated gene or a fragmented one.
Obviously, fragmented genes must be defragmented before  the DNA similarity computation stage (remark that such a defragmentation has already been realized on NCBI website). As the orientation of each fragment is given in the GeneVision output, this defragmentation consists in concatenating all the possible permutations, and only keeping the permutation with the best similarity score to other sequences having the same gene name, if this score is larger than the given threshold.

To put it in a nutshell, the genomes list of gene names are firstly updated in this third approach, following the process detailed in 
Algorithm~\ref{Alg3:thirdM}, while Algorithm~\ref{Alg3:genechk} outlines  
the \emph{geneChk} subroutine. These updated genomes are secondly sent to
Algorithm~\cite{Alkindy2014}, which will produce the desired core genomes, see
Figure~\ref{wholesystem} for an updated pipeline.

%The first method, based  on gene name, translates each genome
% %into a set  of genes without duplicates. 
% The  second method avoid gene
% duplication  through a  defragment  process. In  each iteration,  this
% process  starts by taking  a gene  from gene  list, searches  for gene
% duplication, if a duplication is found, it looks on the orientation of
% the  fragment sequence.   If it  is positive  it appends  directly the
% sequence to  gene files.  Otherwise reverse  complement operations are
% applied  on the sequence,  which is  then also  append to  gene files.
% Finally, a  check for missing start  and stop codons  is performed. At
% the  end  of  the  annotation  process,  all  the  genomes  are  fully
% annotated,  their   genes  are  defragmented,  and   gene  counts  are
% available.
% Remark that there is no gene duplication with gene annotations 
% from DOGMA after applying gene de-fragmentation process. In  
% fact, genome annotation with DOGMA can be the key difference when extracting core genes.
% 
% 
% 
% 
% \subsubsection{Preprocessing}
% 
% In order  to extract  core genomes in  a suitable manner,  the genomic
% data are preprocessed with two methods: on the one hand a method based
% on gene  name and count,  and on  the other hand  a method based  on a
% sequence quality control test.
% 
% In the first method, we extract  a list of genes from each chloroplast
% genome.  Then we store this list of genes in the database under genome
% nam and  genes counts can be  extracted by a  specific length command.
% The \textit{Intersection  Core Matrix}, described  in next subsection,
% is then  computed to  extract the core  genes.  The problem  with this
% method can be stated as follows: how can we ensure that the gene which
% is  predicted in  core genes  is the  same gene  in leaf  genomes? The
% answer  to this problem  is that  if the  sequences of  any gene  in a
% genome annotated  from DOGMA  and NCBI are  similar with respect  to a
% given  threshold,  the method is operational when the sequences are not similar. The problem of attribution of a sequence to a gene in the core genome come to light.
 
% The second method is based on  the underlying idea that it is possible to predict the the best annotated  genome by merging the annotated  genomes from NCBI
% and DOGMA according to a quality test on genes names and sequences. To
% obtain all  quality genes  of each genome,  we consider  the following
% hypothesis: any gene  will appear in the predicted  genome if and only
% if the  annotated genes  in NCBI and  DOGMA pass a  specific threshold
% of  \textit{quality  control test}.   In  fact,  the Needle-man  Wunch
% algorithm  is applied  to compare  both  sequences with  respect to  a
% threshold. If  the alignment  score is above  the threshold,  then the
% gene will be  retained in the predicted genome,  otherwise the gene is
% ignored.   Once    the   prediction   of   all    genomes   is   done,
% the \textit{Intersection Core Matrix} is computed on these new genomes
% to extract core genes, as explained in Algorithm \ref{Alg3:thirdM}.
% 

\begin{algorithm}[H]
\caption{Maximum similarity score between two sequences(geneChk)}
\label{Alg3:genechk}
\begin{algorithmic} 
\REQUIRE $g1,g2 \leftarrow \text{NCBI gene sequence, DOGMA gene sequence}$
\ENSURE $\text{Maximum similarity score}$
\STATE $score1 \leftarrow needle(g1,g2)$
\STATE $score2 \leftarrow needle(g1,Reverse(Complement(g2)))$
\RETURN $max(score1,score2)$
\end{algorithmic}
\end{algorithm}


\begin{algorithm}[H]
\tiny
\caption{Extract new genome based on gene quality test}\label{Alg3:thirdM}
\begin{algorithmic} 
\REQUIRE {$Gname \leftarrow \text{Genome Name}, Threshold \leftarrow 60, RNGenes \leftarrow \text{[ ]}, RDGenes \leftarrow \text{[ ]},PNGenes \leftarrow \text{[ ]}, PDGenes \leftarrow \text{[ ]}$}
\ENSURE $geneList \leftarrow \text{Quality genes}$ 
\FOR{$\text{gene in NCBI genes of Gname}$}
    \IF {$\text{gene in RNGenes}$}
        \STATE $dir(NCBI\_Genes) \leftarrow savePermutation(gene)$
        \STATE $PNGenes \leftarrow gene$
    \ELSE 
    	\STATE $RNGenes \leftarrow gene$
    \ENDIF
\ENDFOR
\FOR{$\text{gene in Dogma genes of Gname}$}
    \IF {$\text{gene in RDGenes}$}
        \STATE $dir(Dogma\_Genes) \leftarrow savePermutation(gene)$
        \STATE $PDGenes \leftarrow gene$
    \ELSE 
    	\STATE $RDGenes \leftarrow gene$
    \ENDIF
\ENDFOR
\STATE $geneList=\text{empty list}$
\STATE $common=set(dir(NCBI\_Genes)) \cap set(dir(Dogma\_Genes))$
\FOR{$\text{gene in common}$}
    \STATE $\text{scores} \leftarrow \text{[ ]}$
    \IF {$\text{gene NOT in PNGenes AND gene NOT in PDGenes}$}
 		\STATE \dots
        \STATE $scores \leftarrow geneChk(g1,g2)$
    \ELSIF {$\text{gene in PNGenes AND NOT gene in PDGenes}$}
        \STATE $PGene \leftarrow loadPermutations('N',gene)$
        \dots
        \FOR {$\text{X in PGene}$}
        	\STATE \dots
	        \STATE $scores \leftarrow geneChk(g1,g2)$
	    \ENDFOR
	 \ELSIF {$\text{gene in PDGenes} AND \text{gene NOT in PNGenes}$} 
	    \STATE $PGene \leftarrow loadPermutations('D',gene)$
	    \STATE \dots
	    \FOR {$\text{X in PGene}$}
	        \STATE \dots
	      	\STATE $scores \leftarrow geneChk(g1,g2)$
	    \ENDFOR
    \ELSIF {$\text{gene in PDGenes} AND \text{gene in PNGenes}$} 
	    \FOR {$\text{X in loadPermutations('N',gene)}$}
	        \FOR {$\text{Y in loadPermutations('D',gene)}$}
	    		\STATE \dots	
	    	  	\STATE $scores \leftarrow geneChk(g1,g2)$
	    	\ENDFOR
	    \ENDFOR
\STATE $score \leftarrow max(scores)$
\IF {$score > Threshold$}
    \STATE $geneList \leftarrow gene$
\ENDIF    
\ENDIF		    
\ENDFOR
\RETURN $geneList$
\end{algorithmic}
\end{algorithm}

%\begin{algorithm}[H]
%\tiny
%\caption{Extract new genome based on gene quality test}\label{Alg3:thirdM}
%\begin{algorithmic} 
%\REQUIRE {$Gname \leftarrow \text{Genome Name}, Threshold \leftarrow 60, RNGenes \leftarrow \text{[ ]}, RDGenes \leftarrow \text{[ ]},PNGenes \leftarrow \text{[ ]}, PDGenes \leftarrow \text{[ ]}$}
%\ENSURE $geneList \leftarrow \text{Quality genes}$ 
%\FOR{$\text{gene in NCBI genes of Gname}$}
%    \IF {$\text{gene in RNGenes}$}
%        \STATE $dir(NCBI\_Genes) \leftarrow savePermutation(gene)$
%        \STATE $PNGenes \leftarrow gene$
%    \ELSE 
%    	\STATE $RNGenes \leftarrow gene$
%    \ENDIF
%\ENDFOR
%\FOR{$\text{gene in Dogma genes of Gname}$}
%    \IF {$\text{gene in RDGenes}$}
%        \STATE $dir(Dogma\_Genes) \leftarrow savePermutation(gene)$
%        \STATE $PDGenes \leftarrow gene$
%    \ELSE 
%    	\STATE $RDGenes \leftarrow gene$
%    \ENDIF
%\ENDFOR
%\STATE $geneList=\text{empty list}$
%\STATE $common=set(dir(NCBI\_Genes)) \cap set(dir(Dogma\_Genes))$
%\FOR{$\text{gene in common}$}
%    \STATE $\text{scores} \leftarrow \text{[ ]}$
%    \IF {$\text{gene NOT in PNGenes AND gene NOT in PDGenes}$}
%        \STATE $g1 \leftarrow open(NCBI\_Genes(gene)).read()$ 	
%        \STATE $g2 \leftarrow open(Dogma\_Genes(gene)).read()$
%        \STATE $score \leftarrow geneChk(g1,g2)$
%        \IF {$score > Threshold$}
%            \STATE $geneList \leftarrow gene$
%        \ENDIF
%    \ELSIF {$\text{gene in PNGenes AND NOT gene in PDGenes}$}
%        \STATE $PGene \leftarrow loadPermutations('N',gene)$
%        \STATE $g2 \leftarrow open(Dogma\_Genes(gene)).read()$
%        \FOR {$\text{X in PGene}$}
%	        \STATE $g1 \leftarrow open(NCBI\_Genes(X)).read()$ 	
%	        \STATE $scores \leftarrow geneChk(g1,g2)$
%	    \ENDFOR
%	    \STATE $score \leftarrow max(scores)$
%	    \IF {$score > Threshold$}
%		    \STATE $geneList \leftarrow gene$
%	    \ENDIF
%	\ELSIF {$\text{gene in PDGenes} AND \text{gene NOT in PNGenes}$} 
%	    \STATE $PGene \leftarrow loadPermutations('D',gene)$
%	    \STATE $g1 \leftarrow open(NCBI\_Genes(gene)).read()$
%	    \FOR {$\text{X in PGene}$}
%	        \STATE $g2 \leftarrow open(Dogma\_Genes(X)).read()$ 	
%	      	\STATE $scores \leftarrow geneChk(g1,g2)$
%	       	\STATE $score \leftarrow max(scores)$
%	    \ENDFOR
%	    \IF {$score > Threshold$}
%	       	\STATE $geneList \leftarrow gene$
%	    \ENDIF
%    \ELSIF {$\text{gene in PDGenes} AND \text{gene in PNGenes}$} 
%	    \FOR {$\text{X in loadPermutations('N',gene)}$}
%	        \FOR {$\text{Y in loadPermutations('D',gene)}$}
%	    		\STATE $g1 \leftarrow open(NCBI\_Genes(X)).read()$
%	            \STATE $g2 \leftarrow open(Dogma\_Genes(Y)).read()$ 	
%	    	  	\STATE $scores \leftarrow geneChk(g1,g2)$
%	    	\ENDFOR
%	    \STATE $score \leftarrow max(scores)$
%	    \IF {$score > Threshold$}
%	       	\STATE $geneList \leftarrow gene$
%	    \ENDIF
%	    \ENDFOR
%	\ENDIF		    
%\ENDFOR
%\RETURN $geneList$
%\end{algorithmic}
%\end{algorithm}




% 
% \textbf{geneChk} is a subroutine used to find the best similarity score between 
% two gene sequences after applying operations like \textit{reverse}, {\it complement}, 
% and {\it reverse complement}. 


  
