BACKGROUND
the first two plant genomes to be sequenced were arabidopsis and rice.
their sequences were obtained by sequencing a minimum tiling path of bacterial artificial chromosome  clones selected from physical maps.
since then, the number of sequenced plant genomes has increased steadily each year, thanks to considerable decrease in costs and increase in throughput of sequencing technologies .
nowadays, most genome assemblies are produced after whole genome shotgun sequencing  using next generation sequencing .
wgs is based on three main steps: i) assembling raw sequence reads into larger sequences called contigs; ii) building bridges between contigs using end-sequenced dna fragments of various lengths  to generate scaffolds; iii) anchoring scaffolds to chromosomes using genetic mapping data to produce pseudo-molecules.
a major challenge is to generate highly contiguous sequence assemblies from short reads in genomes characterized by sequence redundancy, which is a typical situation for plants.
the main source of redundancy is transposable elements  that represent a large part of plant genomes  .
another source of difficulties are paralogous genes  <cit>  resulting from various types of duplications processes including whole genome duplication  that occurred frequently during the evolution of plants  <cit>  or segmental duplication of various sizes.
repeated sequences are often assembled into a single collapsed region during the assembly steps  <cit> .
once created, a collapsed region is linked to multiple other genomic regions leading to conflicts.
automatic assemblers then face two problematic options, either to assemble anyway with a risk to misassemble non-contiguous regions or to prematurely stop the sequence assembly process.
these constraints are exacerbated with short insert-size paired reads since the insert size will not span repeat elements.
conversely, scaffolding with only very large insert size libraries  limits the integration of small scaffolds in the final assembly.
new approaches are continuously being developed to improve genome sequence assemblies.
they include longer read sequencing, high coverage medium and large insert size libraries  <cit> , optical maps , which improve contigs assembly into scaffolds, and genotyping by sequencing , which has been used to assemble scaffolds into pseudo-molecules  <cit> .
in contrast to tremendous advances in high-throughput sequencing, assembling sequences remains a substantial endeavor  <cit> .
several automated programs have been developed to improve draft genome sequence assemblies such as bambus  <cit> , sopra  <cit> , mip  <cit> , sspace  <cit> , opera  <cit> , grass  <cit> , scarpa  <cit> , sspace-longread  <cit> , soap-de-novo <dig>  <cit> , gapfiller  <cit>  and pagit  <cit> .
however, these programs were designed for assembling contigs into scaffolds and/or filling unknown regions, and are running under a compromise between the quantity and quality of the assembly.
this compromise results in a significant proportion of misassembled, un-scaffolded and un-filled regions.
a draft genome sequence assembly of banana , was produced recently using the wgs strategy  <cit> .
the sequence was obtained from a doubled-haploid plant of cv.
pahang and represented a major step forward in understanding the structure and evolution of the banana genome  <cit> .
specific ancestral whole genome duplications were identified within the musa lineage and their impact on gene fractionation and expression patterns was characterized  <cit> .
being the first monocotyledon genome sequence outside the poales, the sequence provided an essential bridge for comparative genome analysis in plants e.g.
.
according to criteria outlined by  <cit> , this genome sequence can be classified as high quality draft.
however, there has been an obvious room for improvement, including the reduction of the number of scaffolds  and the number of scaffolds not anchored to one of the eleven chromosomes .
here we describe a significant improvement of the first musa acuminata draft reference genome sequence and the bioinformatics tools that we developed and used in this work.
the work comprised: i) detection and correction of sequence misassemblies, ii) merging scaffolds, and iii) integration of many previously un-anchored scaffolds to the  <dig> pseudo-molecules.
in addition, conciliation between existing genome annotations was made.
methods
sequence data
the first draft reference sequence of banana   <cit>  was produced from dna of a doubled-haploid plant of cv.
‘pahang’  using reads obtained by  <dig> sequencing , sanger  <dig> kb fosmid paired-reads , sanger bac-end reads  and  <dig> bp pair-end illumina sequences .
in the present work a  <dig> kb mate-pair library of dh-pahang was created and sequenced using illumina hiseq  <dig> to 40x genome coverage.
the reads obtained were trimmed and filtered following three criteria:  trimming of both read ends until base quality is higher or equal to 20;  read trimming at the second unknown base in the sequence; and  read larger or equal to  <dig> bases were conserved.
single molecule mapping
genome map of dh-pahang genome was constructed using bionano irys system .
high molecular weight  dna was prepared according to  <cit> .
briefly, a liquid suspension of intact cell nuclei was prepared by mechanical homogenization of formaldehyde-fixed tissues of unopened  leaves.
the nuclei in the homogenate were stained by dapi , the nuclei in g <dig> phase of cell cycle were purified by flow cytometric sorting and embedded in agarose miniplugs.
hmw dna was then purified and labeled using irysprep reagent kit .
the labelling was done with fluorescent nucleotide analogs at all nt.bspqi nicking endonuclease sites.
single molecules were linearized in nanochannel arrays, imaged.
a total of  <dig>  molecules, with a n <dig> of  <dig> kb, representing a cumulated length of  <dig>  mb with an average label density of  <dig>  labels/ <dig> kb were generated and de novo assembled using a layout-overlap-consensus method.
the de novo map assembly yielded  <dig> mb with a map n <dig> of  <dig> kb.
genetic markers
a total of  <dig> individuals among the  <dig> individuals of a self-progeny of the ‘pahang’ accession  obtained at the cirad research station in guadeloupe were genotyped using the dartseq technology  <cit> .
a total of  <dig>  co-dominant  and  <dig>  dominant markers were generated using a psti-msei enzyme combination.
these markers were used in addition to the  <dig> ssr and  <dig> dart markers previously used to anchor the musa acuminata genome assembly.
out of the  <dig> individuals in the mapping population,  <dig> individuals were genotyped with all types of markers,  <dig> individuals with both dart and dartseq markers,  <dig> individuals with both dartseq and ssr markers and  <dig> individuals with both dart and ssr markers.
the markers were filtered independently for each marker type on the basis of the following criteria: no more than 20 % missing data, no less than 10 % heterozygous or dominant and no less than  <dig>  % homozygous for at least one homozygous state, resulting in  <dig>  markers.
the choice of these relatively non-stringent parameters was motivated by large segregation distortions that were previously observed in chromosome  <dig> and chromosome  <dig> in the segregating population  <cit> .
gene annotation
two gene annotations of the musa acuminata draft genome sequence were available for the initial assembly.
the first corresponded to the annotation published by  <cit> , in addition to approximately  <dig> genes curated by human expertise before  <dig> december  <dig> .
the second one was the ncbi refseq genome annotation released the  <dig> october  <dig>  and generated with the ncbi eukaryotic genome annotation pipeline.
bioinformatics pipeline
an overview of the pipeline used to improve the banana draft genome assembly is shown in fig.
<dig>  it is divided into  <dig> distinct modules corresponding to different and optional operations.
this pipeline exploited several tools that we have developed and which are available under scaffhunter and scaffremodler toolboxes.
the first one exploits genetic mapping data and the second one large insert size paired reads .
they are described in details in the additional file  <dig> fig.
<dig> overview of the pipeline used to improve the musa draft genome sequence.
ellipses correspond to input data and grey ellipses indicate new data acquired for the improvement of the assembly.
boxes corresponds to bioinformatics tools, the ones in blue are new and made available through scaffremodler and scaffhunter toolboxes respectively
module 1: scaffolding of contigs
this module used sspace  <cit>  and exploited large insert size paired reads  to perform a new scaffolding of the existing contigs.
the scaffolding process was divided into as many steps as the number of sequenced libraries with distinct inserts sizes.
the libraries were used by increasing insert size order; scaffolding parameters were optimized for each step.
to prevent accumulation of scaffolding errors, the first library was used with more stringent parameters  than the second and third ones .
for sanger sequence libraries  reads were mapped as single end-reads using bwa  <cit> .
single location reads were used to reconstruct read-pairs that were stored in a tabulated file used by sspace.
module 2: identification and splitting of scaffold/contig misassemblies
this module identified and split misassembled contigs/scaffolds using a combination of gbs genetic mapping data and lpr data.
genetic markers were grouped into linkage groups using joinmap <dig>  software  <cit> .
no marker ordering was performed at this stage.
in parallel, marker sequences were aligned to scaffolds using a consensus of bwa, bowtie <dig>  <cit>  and blast  <cit>  and only single hits markers were conserved.
scaffolds harboring markers attributed to more than one linkage group were identified.
lpr aligning  in these scaffolds were inspected to precisely locate the misassembly boundaries.
the misassembled boundaries were identified based on the absence of overlap of read-pairs in the area and an increased proportion of discordant reads.
misassembled scaffolds were then split.
the complete process and tools used for this module are described in additional file  <dig>
module 3: scaffold fusions/junctions
this module used lpr to identify scaffolds that should be inserted into larger ones  and scaffolds that should be end-joined .
lpr were aligned to the scaffolds using bowtie <dig> in --very-sensitive mode.
only single hit lpr were conserved.
redundant lpr were filtered out using markduplicates tool of picard .
filtered lpr were then used to identify discordant read clusters, which were used to identify potential scaffold fusions and scaffold junctions.
potential scaffold fusions and junctions were then manually validated by inspecting circos  <cit>  picture showing paired reads position in these regions.
fusion and junction performed were validated by aligning lpr along the corrected scaffolds using bowtie <dig>  and mapped reads were inspected to ensure that the newly created junctions are spanned with reads mapped in the correct orientation.
the complete process and tools used for this module are described in additional file  <dig>
module 4: scaffold gap re-estimation
in this module, the size of all remaining gaps  were re-estimated using all paired-reads .
paired reads were aligned against scaffolds using bowtie <dig> in --very-sensitive mode for illumina reads and bwa with mem algorithm for sanger reads.
for each paired read library, gaps were re-estimated so that correctly orientated paired read overlapping a gap have an insert size corresponding to the expected median insert size of the library.
for the  <dig> kb mate-pair library , at least  <dig> pairs were required to re-estimate a gap while for the  <dig> kb and bac-end sanger reads at least  <dig> and  <dig> pairs were required respectively.
the complete process and tools used for this module are described in additional file  <dig>
module 5: super scaffold construction
this module exploits genome map to arrange scaffolds into super scaffolds.
first, the sequence assembly fasta file was converted into the bionano irys map format by running an “in silico digest with the nt.bspqi nicking endonuclease” of the sequence assembly using knickers .
only scaffolds larger than  <dig> kb with more than five sites were used, representing  <dig> scaffolds for a cumulative size of  <dig> mb.
then, using bionano’s proprietary alignment tool refaligner  <cit> , the sequence maps were compared with irys genome maps to find their best alignments; here only sequence maps with more than  <dig> labels  were used for comparison.
the sequence-irys map pairs with significant discordance were flagged and removed, with discordance defined as more than  <dig> consecutive labels not unaligned on both the sequence map and the irys map.
these pairs may represent chimeric assemblies due to sequencing errors or allelic differences.
then, the filtered sequence maps and filtered irys maps were merged with refaligner using a p-value of 10− <dig> based on  <cit>  to create super scaffolds.
this merging process was iterative, and the merge order was based on map similarity.
the iterations stopped when all possible pairs were merged.
a tabulated file locating scaffold sequence into the merged maps was then used to group scaffolds into super scaffolds.
original scaffolds were separated by ns corresponding to their expected distance in the physical map.
module 6: scaffold gap closure
this module exploited paired short insert size reads  to close gaps in scaffold using gapcloser v <dig>  program  <cit> .
at the end of this module, all scaffolds were renamed according to their length.
module 7: scaffold anchoring
this module used genetic markers obtained from a genetic mapping population to group, order and assemble scaffolds into pseudo-molecules.
our approach avoided the step of genetic map construction and a subsequent conciliation between genetic map and scaffolds.
we used blocks of already ordered markers based on their position on scaffolds and first ordered them relative to each other, using upgma-like based approach.
then this first order was improved with permutation testing.
the process can be decomposed into  <dig> steps:marker location on scaffolds using a consensus of bwa, bowtie <dig> and blast,
pairwise linkage lod calculation between markers using joinmap <dig> ,
calculation of a first order using an upgma like approach on mean pairwise linkage lod calculated between scaffolds,
scaffold ordering and orientation optimization by performing scaffold permutations and re-orientations leading to maximization of a score calculated as follows:
score=∑i= <dig> j= <dig> xi<xjn1‐xj‐xinlodij
with n the number of markers in the lg to order, xi and xj are the position of markers i and j in the tested order, and lodij the lod score between markers i and j. to optimize computation time and as order is not tested within scaffolds, i and j are markers from different scaffolds.
scaffold sequences were then assembled into pseudo-molecules.
in addition to a fasta file containing ordered scaffold sequences separated by  <dig> n, an agp file locating scaffolds into pseudo-molecules was generated.
the complete process and tools used for this module are described in additional file  <dig>
module 8: annotation transposition
this module consisted of transposing annotations from the first draft genome sequence to the new assembly.
gene annotations  were transferred to the new assembly using exonerate software  <cit>  with the cdna2genome model and a maximum allowed intron size of  <dig> kb.
exonerate performed genomic searches and spliced alignments in a single run.
using a custom perl script, based on the exonerate output, we transferred the annotation on a new gff <dig> files, and generated a file of sequence identifier equivalence between the two releases.
the script performed some quality checks by comparing protein-coding sequences before and after the transfer as some discrepancies may occur.
in such case, the script used blastp to align genes exons by exons.
since two annotations were available  both annotations were transposed.
an additional consensus annotation was generated using a custom script that selected between the two annotations version genes spanning the same genomic coordinates based on tags enclosed in the gff <dig> files using the intersect function of bedtools  <cit> .
RESULTS
the original banana, musa acuminata, draft reference genome assembly  <cit>  was improved using the approach, tools and datasets as summarized in fig.
<dig>  the improvement was made in  <dig> successive steps.
contig scaffolding
the original  <dig>  contigs published in the first version of the musa acuminata reference genome  <cit>  were re-assembled into scaffolds exploiting paired end data, which were used for original version of the assembly , and new  <dig> kb mate-pair illumina sequences .
contigs were assembled into  <dig>  scaffolds for a cumulated size of  <dig> mb representing 84 % of the estimated size  of the dh-pahang genome .
ninety percent of the assembly was in  <dig> scaffolds and the n <dig> was  <dig>  mb.
gaps  in scaffold represent  <dig>  mb accounting for 11 % of the assembly.table  <dig> statistics on scaffold assemblies
scaffold correction
first, we looked for misassembled scaffolds.
a total of  <dig> scaffolds were identified as containing markers from different linkage groups and thus as potentially containing misassembled regions.
the misassembled regions were confirmed by the presence of discordant  <dig> kb lpr in the region.
the  <dig> misassembled regions identified in these  <dig> scaffolds were then split, resulting in a total of  <dig>  scaffolds.
figure  <dig> shows an example of a misassembled scaffold.
most of the misassembled regions  resulted from scaffolding errors, potentially due to chimeric paired reads or read misalignment.
the remaining misassembled regions  resulted from contig assembly errors.fig.
<dig> example of a clue leading to scaffold splitting.
a genetic markers mapped onto scaffold <dig> belong respectively to linkage-group  <dig>  and linkage-group  <dig>  suggesting a chimeric misassembly.
b circos graphical representation of paired read mapping in the misassembled region.
this representation is drawn using scaffremodler’s tools.
in the inner circle, links between read pairs are drawn with the following color code: grey lines correspond to concordant pairs , orange and red lines correspond to discordant pairs with smaller and greater insert size respectively.
purple lines correspond to pairs showing reverse-reverse orientation, green lines, forward-forward and blue lines correspond to pair with complete reverse orientation relative to the paired library construction.
the second circle represents scaffold in blue with gaps as black regions.
the next circles are scatter plots with warm-cool color code.
the first scatter plot presents the proportion of discordant reads on window size of one third of expected read pair insert size.
the outer circle represents a scatter plot of read coverage on window size of  <dig> bases.
the black arrow points the misassembled region in scaffold <dig> leading to the assembly of two regions that are not linked
second, we looked for potential scaffold fusions and junctions.
based on the analysis of discordant paired-reads from the  <dig> kb lpr with our semi-automated tools, we could perform a total of  <dig> scaffold fusions and  <dig> scaffold junctions, resulting in reduction of scaffold number from  <dig>  to  <dig> .
figure  <dig> shows an example of clue leading to scaffold <dig> fusion into scaffold <dig>  additional file 1: figure s <dig> shows the mapping of reads on the two borders of scaffold <dig> after fusion into scaffold <dig>  both right and left borders displayed overlapping reads in the correct orientation .fig.
<dig> example of a clue leading to scaffold fusion.
a graphical representation of paired read leading to the identification of fusion of scaffold <dig> into scaffold <dig>  this representation is drawn using scaffremodler’s tools.
in the inner circle, links between read pairs are drawn with the color code described in fig.
2: grey for concordant pairs; red and orange for discordant in size; purple, green and blue for orientation discordant pairs.
the second circle represents scaffold in blue with gaps as black regions.
the next represents the proportion of discordant reads and the last circle represents read coverage as in fig.
<dig>  red and blue beams linking scaffold <dig> and scaffold <dig> allowed identifying scaffold fusion schematized in .
inserting scaffold <dig> into scaffold <dig> will correct the discordant red links and correct the orientation of discordant blue links
at this stage the size of gaps  within the new  <dig>  scaffolds was re-estimated using the paired reads libraries sequentially resulting in  <dig> mb for  <dig>  % of the assembly .
the cumulative size of the new  <dig>  scaffolds after gap re-estimation was of  <dig> mb.
ninety percent of the assembly was in  <dig> scaffolds and the n <dig> was  <dig>  mb.
finally, bionano irys genome map of dh-pahang was used to order and orient scaffolds into super scaffolds.
this step allowed merging of  <dig> scaffolds into  <dig> super-scaffolds.
a total of  <dig>  mb of gap regions were added during super scaffold construction .
finally, 90 % of the assembly was in  <dig> scaffolds and the n <dig> was  <dig>  mb with  <dig> scaffolds.
gaps in scaffolds represented  <dig>  mb for  <dig>  % of the assembly.
gap closure
gaps within the  <dig>  scaffolds were then tentatively filled with the gapcloser program using the  <dig> bp pair-end illumina sequencing libraries , generated to correct the first version of the banana musa acuminata reference genome.
of the total of  <dig>  gap regions,  <dig>  were closed.
final assembly
the final assembly  consisted of  <dig>  scaffolds and showed a cumulative size of  <dig>  mb corresponding to 86 % of the estimated size of the dh-pahang genome.
ninety percent of the assembly was in  <dig> scaffolds and the n <dig> was  <dig>  mb.
gaps in scaffolds represent only  <dig>  mb .
twelve of these scaffolds were identified as mitochondrial dna  using blast  of mitochondrial coding sequences of phoenix dactylifera .
the twelve mitochondrial scaffolds were removed from the final nuclear assembly.
in order to validate the improvements made, the proportion of mapped  <dig> kb mate pair discordant reads  for each scaffold assembly versions was calculated.
over the  <dig>  million non-redundant and single mapped pairs,  <dig>  million  mapped discordantly on the first version.
over the  <dig>  million non-redundant and single mapped pairs  <dig>  million  mapped discordantly for the new assembly before gap closure.
over the  <dig>  million non-redundant and single mapped pairs  <dig>  million  mapped discordantly for the new assembly after gap closure.
musa scaffold anchoring
genetic markers were then used to assemble scaffolds into pseudo-molecules.
of the  <dig>  selected genetic markers,  <dig>  that mapped to a unique position were grouped into  <dig> linkage groups.
a total of  <dig> markers were discarded since they created local discrepancies in scaffolds, clearly attributed to a linkage group based on the majority of the markers.
markers located on small scaffolds for which no linkage group majority could be found were also discarded.
the remaining  <dig>  markers allowed to order and orient  <dig> scaffolds into the  <dig> pseudo-molecules , with an average of  <dig>  markers per  <dig> kb .fig.
<dig> representation of the new version of eleven pseudo-molecules of musa acuminata.
black and white boxes correspond to oriented and unoriented scaffolds, respectively.
genetic marker, gene and unknown sequence  density are represented in grey, blue and green respectively based on a windows size of  <dig> kb.
the recombination rate  has been calculated on  <dig> individuals on corrected genetic markers and a sliding window of  <dig> kb
finally, a total of  <dig> mb of genome sequence was anchored, representing  <dig>  % of the nuclear genome assembly  and including all scaffolds larger than  <dig> mb.
each pseudo-molecule comprised from  <dig> to  <dig> scaffolds and n <dig> in pseudo-molecules varied between  <dig>  mb to  <dig>  mb.
the mean n  proportion varied from  <dig>  to  <dig>  % in pseudo-molecules and was of  <dig>  % in unanchored scaffolds .
marker linkage in ordered scaffolds can be visualized for each chromosome in additional file 1: figure s <dig> table  <dig> statistics on musa acuminata pseudo-molecule assembly between the first and the new version
ascaffold number
in comparison to the first pseudo-molecule assembly version, we corrected the position of only a few large regions from one pseudo-molecule to another .
one major change concerned a region that was previously anchored to chromosome  <dig> and that is now assigned to chromosome  <dig>  these regions of chromosomes  <dig> and  <dig> displayed marked segregation distortions that created pseudo-linkages  <cit>  and hampered the anchoring of the first draft assembly that was based on much lower number of genetic markers.
apart from this large change in the assembly, many small modifications were made, representing either anchoring small scaffolds previously unanchored, or small scaffolds reordering.
most of these changes concerned peri-centromeric regions.fig.
<dig> dot plot comparison of gene order between the initial and the new version of musa acuminata genome sequence assembly.
a dot represents the position of a gene in the two assembly versions with the initial assembly on x axis and the new one on the y axis.
ruptures in the diagonal indicate differences of gene order.
red circles indicate the main differences and green circles indicate the variations resulting from the approximate scaffold order in the peri-centromeric regions.
for instance, the version  <dig> of the assembly corrects a significant error between the chromosome  <dig> and 4
annotation transfer
two independent annotations of the initial version of the banana genome assembly were available and both were transferred to the new assembly.
the m. acuminata transcripts from the first annotation published  <cit>  in addition to several manually curated gene annotation were transferred to the new assembly version.
of the  <dig>  predicted genes,  <dig>   genes were transferred to the new assembly version .
of the total number of transferred genes,  <dig>  were located in unanchored scaffolds compared to  <dig>  genes  in the first version.
ninety-six genes were transferred onto the mitochondrial scaffolds.
the same transfer was performed for the ncbi refseq genome annotation.
a total of  <dig>   genes of the  <dig>  predicted genes were transferred to the new assembly version .table  <dig> statistics on annotation transfer between the first release of the assembly and the new release
ancbi refseq genome annotation released the  <dig> october  <dig> and generated with the ncbi eukaryotic genome annotation pipeline
bbanana genome hub  annotation performed by  <cit> , in addition to manually curated genes performed before  <dig> december  <dig> available in the banana genome hub
cincluding ‘n’ separating scaffolds
based on the analysis of several manually curated genes, the ncbi refseq genome annotation proved to be generally of better quality than the first published annotation in particular because the first annotation over predicted introns.
in addition, the ncbi refseq genome annotation integrated rnaseq data and predicted alternative transcripts.
we thus created a consensus annotation that combined all the manually curated genes, the ncbi refseq annotation and the predicted genes from the first annotation that were missed by the refseq annotation pipeline.
using jbrowse in the banana genome hub, these three gene annotations can be visualized as separate tracks.
note that since, we did not perform a new annotation but an annotation transfer, gene fragmentation due to contigs miss-junctions still remains in the new annotated assembly version even if the new assembly version corrected such gene fragmentation.
finally, the consensus annotation contains  <dig>  predicted genes with  <dig>   located in chromosomes,  <dig>  located in unanchored scaffolds and  <dig>  located in identified mitochondrial scaffolds .
to avoid any confusion, we modified the nomenclature of locus tags.
for example, gsmua_achr5t02570_ <dig> in version  <dig> becomes ma05_t <dig>  in version  <dig>
discussion
during the course of this work we succeeded in significantly improving the initial musa nuclear draft genome assembly by reducing the scaffold number by 80 % , doubling the n <dig> value  and increasing the proportion of assembly anchored to the  <dig> musa chromosomes by 20 %  that now include  <dig>  % of genes.
the decrease of discordant  <dig> kb read-pairs mapping proportion of 40 % between initial and new version of the assembly support the quality of the changes that were made.
the addition of the  <dig> kb mate-pair illumina library in the scaffolding process decreased scaffold number by 70 %  and raised n <dig> from  <dig>  mb to  <dig>  mb.
these results highlighted the importance of medium insert size library during the scaffolding process.
interestingly, the scaffold fusion/junction that we performed decreased further the scaffold number by 30 %  and significantly impacted the n <dig> value which nearly doubled.
these results highlight the utility and power of the semi-automated tools we have developed.
apart from verifying the newly established scaffolds, the use of bionano irys genome maps permitted a few additional scaffold junctions.
these maps would have had a bigger impact if they were available earlier during the process  <cit> .
the gap filling step allowed an important reduction of gap regions in the final assembly .
the reduction of discordant  <dig> kb read pairs proportion between the assembly before and after gap filling highlighted the quality of gap closure step performed.
the cumulative size of the new assembly is reduced by  <dig>  mb in comparison with the first genome assembly  <cit> .
this reduction is mainly due to the insertion of small scaffolds into previous gaps of larger scaffolds.
the total size of the assembly, lower than expected, can be explained at least in part by difficulties in correctly assembling the repeated fraction of the genome .
these repeat-rich sequences are often collapsed into single regions, resulting in a reduced size for the total assembly  <cit> .
for example,  <dig>  mb rdna have been found in the unassembled reads of dh-pahang  <cit> .
saturation of genetic map with dartseq markers increased the proportion of anchored assembly from  <dig> to  <dig>  % and anchored genes from  <dig> to  <dig>  %.
for scaffold anchoring, the classical approach is to construct a genetic map and to anchor the scaffold assembly onto this genetic map to construct a pseudo-molecule.
genotyping errors that are frequent in gbs data can lead to marker miss-ordering in genetic map and to conflict between markers order in genetic map and in scaffolds, when performing the scaffold anchoring.
to avoid the tedious step of conciliation between genetic map and scaffolds, we developed a method that takes the advantage of markers already ordered into blocks corresponding to scaffolds.
in this context, genotyping error impact is lowered as markers are already partially ordered.
the newly anchored regions belong essentially to peri-centromeric regions.
however because the proportion of repeated sequence is high in these regions, the marker density is lower  and the recombination rate is generally very low  .
consequently the scaffold order and orientation in these regions remains tentative.
CONCLUSIONS
the significant improvements made on the banana reference genome sequence will have important impact on the quality of future genetic and comparative genomic analysis.
the bioinformatics methods and tools described in this work can be useful to improve draft genome assemblies in other plant species.
the pipeline comprises independent modules adaptable to various datatypes.
it can be used to improve existing assemblies or in combination with existing automated programs during de novo assembly.
the improved version of the musa acuminata genome assembly is accessible and can be downloaded in the new version of the banana genome hub at http://banana-genome.cirad.fr/  <cit> .
tools are available in command line version on github .
most of the options  are also available on the south green galaxy platform under scaffhunter and scaffremodler toolboxes .
availability of supporting data
datasets  are available through the banana genome hub  and the  <dig> kb library is deposited on the ena read archive .
additional file
additional file 1: detailed description of tools and processes used to improve the musa acuminata reference sequence and additional figures.
abbreviations
bacbacterial artificial chromosome
gbsgenotyping by sequencing
hmwhigh molecular weight
lprlarge insert size paired reads
ngsnext generation sequencing
tetransposable elements
wgdwhole genome duplication
wgswhole genome shotgun sequencing
competing interests
the authors declare that they have no competing interests.
authors’ contributions
gm, fcb, adh: conceived and designed the study and wrote the manuscript.
gm: developed the bioinformatic programs and performed the analysis.
ak, aa, jd, ah: produced the sequencing data and the genome irys map.
gd, mr, ac, ak, jma, ah, fc: contributed to the analysis and edited the manuscript.
adh: coordinated the study.
all authors read and approved the final manuscript.
the authors thank the diversity arrays technology pty ltd for dartseq genotyping, jan vrána and hana Šimková for preparation of hmw dna, and cgiar research program on roots, tubers and bananas  for financial support for sequencing data acquisition.
we also thank the south green bioinformatics platform  for providing us with computational resources.
we thank christophe jenny for providing the pahang segregating population from the cirad research station in guadeloupe, french west indies.

 the first two plant genomes to be sequenced were arabidopsis and rice. their sequences were obtained by sequencing a minimum tiling path of bacterial artificial chromosome clones selected from physical maps. since then, the number of sequenced plant genomes has increased steadily each year, thanks to considerable decrease in costs and increase in throughput of sequencing technologies. nowadays, most genome assemblies are produced after whole genome shotgun sequencing using next generation sequencing. wgs is based on three main steps: i) assembling raw sequence reads into larger sequences called contigs; ii) building bridges between contigs using end-sequenced dna fragments of various lengths to generate scaffolds; iii) anchoring scaffolds to chromosomes using genetic mapping data to produce pseudo-molecules. a major challenge is to generate highly contiguous sequence assemblies from short reads in genomes characterized by sequence redundancy, which is a typical situation for plants. the main source of redundancy is transposable elements that represent a large part of plant genomes. another source of difficulties are paralogous genes resulting from various types of duplications processes including whole genome duplication that occurred frequently during the evolution of plants or segmental duplication of various sizes. repeated sequences are often assembled into a single collapsed region during the assembly steps. once created, a collapsed region is linked to multiple other genomic regions leading to conflicts. automatic assemblers then face two problematic options, either to assemble anyway with a risk to misassemble non-contiguous regions or to prematurely stop the sequence assembly process. these constraints are exacerbated with short insert-size paired reads since the insert size will not span repeat elements. conversely, scaffolding with only very large insert size libraries limits the integration of small scaffolds in the final assembly. new approaches are continuously being developed to improve genome sequence assemblies. they include longer read sequencing, high coverage medium and large insert size libraries, optical maps, which improve contigs assembly into scaffolds, and genotyping by sequencing, which has been used to assemble scaffolds into pseudo-molecules. in contrast to tremendous advances in high-throughput sequencing, assembling sequences remains a substantial endeavor. several automated programs have been developed to improve draft genome sequence assemblies such as bambus, sopra, mip, sspace, opera, grass, scarpa, sspace-longread, soap-de-novo, gapfiller and pagit. however, these programs were designed for assembling contigs into scaffolds and/or filling unknown regions, and are running under a compromise between the quantity and quality of the assembly. this compromise results in a significant proportion of misassembled, un-scaffolded and un-filled regions. a draft genome sequence assembly of banana, was produced recently using the wgs strategy. the sequence was obtained from a doubled-haploid plant of cv. pahang and represented a major step forward in understanding the structure and evolution of the banana genome. specific ancestral whole genome duplications were identified within the musa lineage and their impact on gene fractionation and expression patterns was characterized. being the first monocotyledon genome sequence outside the poales, the sequence provided an essential bridge for comparative genome analysis in plants e.g.. according to criteria outlined by, this genome sequence can be classified as high quality draft. however, there has been an obvious room for improvement, including the reduction of the number of scaffolds and the number of scaffolds not anchored to one of the eleven chromosomes. here we describe a significant improvement of the first musa acuminata draft reference genome sequence and the bioinformatics tools that we developed and used in this work. the work comprised: i) detection and correction of sequence misassemblies, ii) merging scaffolds, and iii) integration of many previously un-anchored scaffolds to the pseudo-molecules. in addition, conciliation between existing genome annotations was made. methods sequence data the first draft reference sequence of banana was produced from dna of a doubled-haploid plant of cv. ‘pahang’ using reads obtained by sequencing, sanger kb fosmid paired-reads, sanger bac-end reads and bp pair-end illumina sequences. in the present work a kb mate-pair library of dh-pahang was created and sequenced using illumina hiseq to 40x genome coverage. the reads obtained were trimmed and filtered following three criteria: trimming of both read ends until base quality is higher or equal to 20; read trimming at the second unknown base in the sequence; and read larger or equal to bases were conserved. single molecule mapping genome map of dh-pahang genome was constructed using bionano irys system. high molecular weight dna was prepared according to. briefly, a liquid suspension of intact cell nuclei was prepared by mechanical homogenization of formaldehyde-fixed tissues of unopened leaves. the nuclei in the homogenate were stained by dapi, the nuclei in g phase of cell cycle were purified by flow cytometric sorting and embedded in agarose miniplugs. hmw dna was then purified and labeled using irysprep reagent kit. the labelling was done with fluorescent nucleotide analogs at all nt.bspqi nicking endonuclease sites. single molecules were linearized in nanochannel arrays, imaged. a total of molecules, with a n of kb, representing a cumulated length of mb with an average label density of labels/ kb were generated and de novo assembled using a layout-overlap-consensus method. the de novo map assembly yielded mb with a map n of kb. genetic markers a total of individuals among the individuals of a self-progeny of the ‘pahang’ accession obtained at the cirad research station in guadeloupe were genotyped using the dartseq technology. a total of co-dominant and dominant markers were generated using a psti-msei enzyme combination. these markers were used in addition to the ssr and dart markers previously used to anchor the musa acuminata genome assembly. out of the individuals in the mapping population, individuals were genotyped with all types of markers, individuals with both dart and dartseq markers, individuals with both dartseq and ssr markers and individuals with both dart and ssr markers. the markers were filtered independently for each marker type on the basis of the following criteria: no more than 20 % missing data, no less than 10 % heterozygous or dominant and no less than  % homozygous for at least one homozygous state, resulting in markers. the choice of these relatively non-stringent parameters was motivated by large segregation distortions that were previously observed in chromosome and chromosome in the segregating population. gene annotation two gene annotations of the musa acuminata draft genome sequence were available for the initial assembly. the first corresponded to the annotation published by, in addition to approximately genes curated by human expertise before december. the second one was the ncbi refseq genome annotation released the october and generated with the ncbi eukaryotic genome annotation pipeline. bioinformatics pipeline an overview of the pipeline used to improve the banana draft genome assembly is shown in fig.  it is divided into distinct modules corresponding to different and optional operations. this pipeline exploited several tools that we have developed and which are available under scaffhunter and scaffremodler toolboxes. the first one exploits genetic mapping data and the second one large insert size paired reads. they are described in details in the additional file fig. overview of the pipeline used to improve the musa draft genome sequence. ellipses correspond to input data and grey ellipses indicate new data acquired for the improvement of the assembly. boxes corresponds to bioinformatics tools, the ones in blue are new and made available through scaffremodler and scaffhunter toolboxes respectively module 1: scaffolding of contigs this module used sspace and exploited large insert size paired reads to perform a new scaffolding of the existing contigs. the scaffolding process was divided into as many steps as the number of sequenced libraries with distinct inserts sizes. the libraries were used by increasing insert size order; scaffolding parameters were optimized for each step. to prevent accumulation of scaffolding errors, the first library was used with more stringent parameters than the second and third ones. for sanger sequence libraries reads were mapped as single end-reads using bwa. single location reads were used to reconstruct read-pairs that were stored in a tabulated file used by sspace. module 2: identification and splitting of scaffold/contig misassemblies this module identified and split misassembled contigs/scaffolds using a combination of gbs genetic mapping data and lpr data. genetic markers were grouped into linkage groups using joinmap software. no marker ordering was performed at this stage. in parallel, marker sequences were aligned to scaffolds using a consensus of bwa, bowtie and blast and only single hits markers were conserved. scaffolds harboring markers attributed to more than one linkage group were identified. lpr aligning in these scaffolds were inspected to precisely locate the misassembly boundaries. the misassembled boundaries were identified based on the absence of overlap of read-pairs in the area and an increased proportion of discordant reads. misassembled scaffolds were then split. the complete process and tools used for this module are described in additional file module 3: scaffold fusions/junctions this module used lpr to identify scaffolds that should be inserted into larger ones and scaffolds that should be end-joined. lpr were aligned to the scaffolds using bowtie in --very-sensitive mode. only single hit lpr were conserved. redundant lpr were filtered out using markduplicates tool of picard. filtered lpr were then used to identify discordant read clusters, which were used to identify potential scaffold fusions and scaffold junctions. potential scaffold fusions and junctions were then manually validated by inspecting circos picture showing paired reads position in these regions. fusion and junction performed were validated by aligning lpr along the corrected scaffolds using bowtie and mapped reads were inspected to ensure that the newly created junctions are spanned with reads mapped in the correct orientation. the complete process and tools used for this module are described in additional file module 4: scaffold gap re-estimation in this module, the size of all remaining gaps were re-estimated using all paired-reads. paired reads were aligned against scaffolds using bowtie in --very-sensitive mode for illumina reads and bwa with mem algorithm for sanger reads. for each paired read library, gaps were re-estimated so that correctly orientated paired read overlapping a gap have an insert size corresponding to the expected median insert size of the library. for the kb mate-pair library, at least pairs were required to re-estimate a gap while for the kb and bac-end sanger reads at least and pairs were required respectively. the complete process and tools used for this module are described in additional file module 5: super scaffold construction this module exploits genome map to arrange scaffolds into super scaffolds. first, the sequence assembly fasta file was converted into the bionano irys map format by running an “in silico digest with the nt.bspqi nicking endonuclease” of the sequence assembly using knickers. only scaffolds larger than kb with more than five sites were used, representing scaffolds for a cumulative size of mb. then, using bionano’s proprietary alignment tool refaligner, the sequence maps were compared with irys genome maps to find their best alignments; here only sequence maps with more than labels were used for comparison. the sequence-irys map pairs with significant discordance were flagged and removed, with discordance defined as more than consecutive labels not unaligned on both the sequence map and the irys map. these pairs may represent chimeric assemblies due to sequencing errors or allelic differences. then, the filtered sequence maps and filtered irys maps were merged with refaligner using a p-value of 10− based on to create super scaffolds. this merging process was iterative, and the merge order was based on map similarity. the iterations stopped when all possible pairs were merged. a tabulated file locating scaffold sequence into the merged maps was then used to group scaffolds into super scaffolds. original scaffolds were separated by ns corresponding to their expected distance in the physical map. module 6: scaffold gap closure this module exploited paired short insert size reads to close gaps in scaffold using gapcloser v program. at the end of this module, all scaffolds were renamed according to their length. module 7: scaffold anchoring this module used genetic markers obtained from a genetic mapping population to group, order and assemble scaffolds into pseudo-molecules. our approach avoided the step of genetic map construction and a subsequent conciliation between genetic map and scaffolds. we used blocks of already ordered markers based on their position on scaffolds and first ordered them relative to each other, using upgma-like based approach. then this first order was improved with permutation testing. the process can be decomposed into steps:marker location on scaffolds using a consensus of bwa, bowtie and blast, pairwise linkage lod calculation between markers using joinmap, calculation of a first order using an upgma like approach on mean pairwise linkage lod calculated between scaffolds, scaffold ordering and orientation optimization by performing scaffold permutations and re-orientations leading to maximization of a score calculated as follows: score=∑i= j= xi<xjn1‐xj‐xinlodij with n the number of markers in the lg to order, xi and xj are the position of markers i and j in the tested order, and lodij the lod score between markers i and j. to optimize computation time and as order is not tested within scaffolds, i and j are markers from different scaffolds. scaffold sequences were then assembled into pseudo-molecules. in addition to a fasta file containing ordered scaffold sequences separated by n, an agp file locating scaffolds into pseudo-molecules was generated. the complete process and tools used for this module are described in additional file module 8: annotation transposition this module consisted of transposing annotations from the first draft genome sequence to the new assembly. gene annotations were transferred to the new assembly using exonerate software with the cdna2genome model and a maximum allowed intron size of kb. exonerate performed genomic searches and spliced alignments in a single run. using a custom perl script, based on the exonerate output, we transferred the annotation on a new gff files, and generated a file of sequence identifier equivalence between the two releases. the script performed some quality checks by comparing protein-coding sequences before and after the transfer as some discrepancies may occur. in such case, the script used blastp to align genes exons by exons. since two annotations were available both annotations were transposed. an additional consensus annotation was generated using a custom script that selected between the two annotations version genes spanning the same genomic coordinates based on tags enclosed in the gff files using the intersect function of bedtools.  the original banana, musa acuminata, draft reference genome assembly was improved using the approach, tools and datasets as summarized in fig.  the improvement was made in successive steps. contig scaffolding the original contigs published in the first version of the musa acuminata reference genome were re-assembled into scaffolds exploiting paired end data, which were used for original version of the assembly, and new kb mate-pair illumina sequences. contigs were assembled into scaffolds for a cumulated size of mb representing 84 % of the estimated size of the dh-pahang genome. ninety percent of the assembly was in scaffolds and the n was mb. gaps in scaffold represent mb accounting for 11 % of the assembly.table statistics on scaffold assemblies scaffold correction first, we looked for misassembled scaffolds. a total of scaffolds were identified as containing markers from different linkage groups and thus as potentially containing misassembled regions. the misassembled regions were confirmed by the presence of discordant kb lpr in the region. the misassembled regions identified in these scaffolds were then split, resulting in a total of scaffolds. figure  shows an example of a misassembled scaffold. most of the misassembled regions resulted from scaffolding errors, potentially due to chimeric paired reads or read misalignment. the remaining misassembled regions resulted from contig assembly errors.fig. example of a clue leading to scaffold splitting. a genetic markers mapped onto scaffold belong respectively to linkage-group and linkage-group suggesting a chimeric misassembly. b circos graphical representation of paired read mapping in the misassembled region. this representation is drawn using scaffremodler’s tools. in the inner circle, links between read pairs are drawn with the following color code: grey lines correspond to concordant pairs, orange and red lines correspond to discordant pairs with smaller and greater insert size respectively. purple lines correspond to pairs showing reverse-reverse orientation, green lines, forward-forward and blue lines correspond to pair with complete reverse orientation relative to the paired library construction. the second circle represents scaffold in blue with gaps as black regions. the next circles are scatter plots with warm-cool color code. the first scatter plot presents the proportion of discordant reads on window size of one third of expected read pair insert size. the outer circle represents a scatter plot of read coverage on window size of bases. the black arrow points the misassembled region in scaffold leading to the assembly of two regions that are not linked second, we looked for potential scaffold fusions and junctions. based on the analysis of discordant paired-reads from the kb lpr with our semi-automated tools, we could perform a total of scaffold fusions and scaffold junctions, resulting in reduction of scaffold number from to. figure  shows an example of clue leading to scaffold fusion into scaffold additional file 1: figure s shows the mapping of reads on the two borders of scaffold after fusion into scaffold both right and left borders displayed overlapping reads in the correct orientation.fig. example of a clue leading to scaffold fusion. a graphical representation of paired read leading to the identification of fusion of scaffold into scaffold this representation is drawn using scaffremodler’s tools. in the inner circle, links between read pairs are drawn with the color code described in fig. 2: grey for concordant pairs; red and orange for discordant in size; purple, green and blue for orientation discordant pairs. the second circle represents scaffold in blue with gaps as black regions. the next represents the proportion of discordant reads and the last circle represents read coverage as in fig.  red and blue beams linking scaffold and scaffold allowed identifying scaffold fusion schematized in. inserting scaffold into scaffold will correct the discordant red links and correct the orientation of discordant blue links at this stage the size of gaps within the new scaffolds was re-estimated using the paired reads libraries sequentially resulting in mb for  % of the assembly. the cumulative size of the new scaffolds after gap re-estimation was of mb. ninety percent of the assembly was in scaffolds and the n was mb. finally, bionano irys genome map of dh-pahang was used to order and orient scaffolds into super scaffolds. this step allowed merging of scaffolds into super-scaffolds. a total of mb of gap regions were added during super scaffold construction. finally, 90 % of the assembly was in scaffolds and the n was mb with scaffolds. gaps in scaffolds represented mb for  % of the assembly. gap closure gaps within the scaffolds were then tentatively filled with the gapcloser program using the bp pair-end illumina sequencing libraries, generated to correct the first version of the banana musa acuminata reference genome. of the total of gap regions, were closed. final assembly the final assembly consisted of scaffolds and showed a cumulative size of mb corresponding to 86 % of the estimated size of the dh-pahang genome. ninety percent of the assembly was in scaffolds and the n was mb. gaps in scaffolds represent only mb. twelve of these scaffolds were identified as mitochondrial dna using blast of mitochondrial coding sequences of phoenix dactylifera. the twelve mitochondrial scaffolds were removed from the final nuclear assembly. in order to validate the improvements made, the proportion of mapped kb mate pair discordant reads for each scaffold assembly versions was calculated. over the million non-redundant and single mapped pairs, million mapped discordantly on the first version. over the million non-redundant and single mapped pairs million mapped discordantly for the new assembly before gap closure. over the million non-redundant and single mapped pairs million mapped discordantly for the new assembly after gap closure. musa scaffold anchoring genetic markers were then used to assemble scaffolds into pseudo-molecules. of the selected genetic markers, that mapped to a unique position were grouped into linkage groups. a total of markers were discarded since they created local discrepancies in scaffolds, clearly attributed to a linkage group based on the majority of the markers. markers located on small scaffolds for which no linkage group majority could be found were also discarded. the remaining markers allowed to order and orient scaffolds into the pseudo-molecules, with an average of markers per kb.fig. representation of the new version of eleven pseudo-molecules of musa acuminata. black and white boxes correspond to oriented and unoriented scaffolds, respectively. genetic marker, gene and unknown sequence density are represented in grey, blue and green respectively based on a windows size of kb. the recombination rate has been calculated on individuals on corrected genetic markers and a sliding window of kb finally, a total of mb of genome sequence was anchored, representing  % of the nuclear genome assembly and including all scaffolds larger than mb. each pseudo-molecule comprised from to scaffolds and n in pseudo-molecules varied between mb to mb. the mean n proportion varied from to  % in pseudo-molecules and was of  % in unanchored scaffolds. marker linkage in ordered scaffolds can be visualized for each chromosome in additional file 1: figure s table statistics on musa acuminata pseudo-molecule assembly between the first and the new version ascaffold number in comparison to the first pseudo-molecule assembly version, we corrected the position of only a few large regions from one pseudo-molecule to another. one major change concerned a region that was previously anchored to chromosome and that is now assigned to chromosome these regions of chromosomes and displayed marked segregation distortions that created pseudo-linkages and hampered the anchoring of the first draft assembly that was based on much lower number of genetic markers. apart from this large change in the assembly, many small modifications were made, representing either anchoring small scaffolds previously unanchored, or small scaffolds reordering. most of these changes concerned peri-centromeric regions.fig. dot plot comparison of gene order between the initial and the new version of musa acuminata genome sequence assembly. a dot represents the position of a gene in the two assembly versions with the initial assembly on x axis and the new one on the y axis. ruptures in the diagonal indicate differences of gene order. red circles indicate the main differences and green circles indicate the variations resulting from the approximate scaffold order in the peri-centromeric regions. for instance, the version of the assembly corrects a significant error between the chromosome and 4 annotation transfer two independent annotations of the initial version of the banana genome assembly were available and both were transferred to the new assembly. the m. acuminata transcripts from the first annotation published in addition to several manually curated gene annotation were transferred to the new assembly version. of the predicted genes, genes were transferred to the new assembly version. of the total number of transferred genes, were located in unanchored scaffolds compared to genes in the first version. ninety-six genes were transferred onto the mitochondrial scaffolds. the same transfer was performed for the ncbi refseq genome annotation. a total of genes of the predicted genes were transferred to the new assembly version.table statistics on annotation transfer between the first release of the assembly and the new release ancbi refseq genome annotation released the october and generated with the ncbi eukaryotic genome annotation pipeline bbanana genome hub annotation performed by, in addition to manually curated genes performed before december available in the banana genome hub cincluding ‘n’ separating scaffolds based on the analysis of several manually curated genes, the ncbi refseq genome annotation proved to be generally of better quality than the first published annotation in particular because the first annotation over predicted introns. in addition, the ncbi refseq genome annotation integrated rnaseq data and predicted alternative transcripts. we thus created a consensus annotation that combined all the manually curated genes, the ncbi refseq annotation and the predicted genes from the first annotation that were missed by the refseq annotation pipeline. using jbrowse in the banana genome hub, these three gene annotations can be visualized as separate tracks. note that since, we did not perform a new annotation but an annotation transfer, gene fragmentation due to contigs miss-junctions still remains in the new annotated assembly version even if the new assembly version corrected such gene fragmentation. finally, the consensus annotation contains predicted genes with located in chromosomes, located in unanchored scaffolds and located in identified mitochondrial scaffolds. to avoid any confusion, we modified the nomenclature of locus tags. for example, gsmua_achr5t02570_ in version becomes ma05_t in version discussion during the course of this work we succeeded in significantly improving the initial musa nuclear draft genome assembly by reducing the scaffold number by 80 %, doubling the n value and increasing the proportion of assembly anchored to the musa chromosomes by 20 % that now include  % of genes. the decrease of discordant kb read-pairs mapping proportion of 40 % between initial and new version of the assembly support the quality of the changes that were made. the addition of the kb mate-pair illumina library in the scaffolding process decreased scaffold number by 70 % and raised n from mb to mb. these results highlighted the importance of medium insert size library during the scaffolding process. interestingly, the scaffold fusion/junction that we performed decreased further the scaffold number by 30 % and significantly impacted the n value which nearly doubled. these results highlight the utility and power of the semi-automated tools we have developed. apart from verifying the newly established scaffolds, the use of bionano irys genome maps permitted a few additional scaffold junctions. these maps would have had a bigger impact if they were available earlier during the process. the gap filling step allowed an important reduction of gap regions in the final assembly. the reduction of discordant kb read pairs proportion between the assembly before and after gap filling highlighted the quality of gap closure step performed. the cumulative size of the new assembly is reduced by mb in comparison with the first genome assembly. this reduction is mainly due to the insertion of small scaffolds into previous gaps of larger scaffolds. the total size of the assembly, lower than expected, can be explained at least in part by difficulties in correctly assembling the repeated fraction of the genome. these repeat-rich sequences are often collapsed into single regions, resulting in a reduced size for the total assembly. for example, mb rdna have been found in the unassembled reads of dh-pahang. saturation of genetic map with dartseq markers increased the proportion of anchored assembly from to  % and anchored genes from to  %. for scaffold anchoring, the classical approach is to construct a genetic map and to anchor the scaffold assembly onto this genetic map to construct a pseudo-molecule. genotyping errors that are frequent in gbs data can lead to marker miss-ordering in genetic map and to conflict between markers order in genetic map and in scaffolds, when performing the scaffold anchoring. to avoid the tedious step of conciliation between genetic map and scaffolds, we developed a method that takes the advantage of markers already ordered into blocks corresponding to scaffolds. in this context, genotyping error impact is lowered as markers are already partially ordered. the newly anchored regions belong essentially to peri-centromeric regions. however because the proportion of repeated sequence is high in these regions, the marker density is lower and the recombination rate is generally very low. consequently the scaffold order and orientation in these regions remains tentative.  the significant improvements made on the banana reference genome sequence will have important impact on the quality of future genetic and comparative genomic analysis. the bioinformatics methods and tools described in this work can be useful to improve draft genome assemblies in other plant species. the pipeline comprises independent modules adaptable to various datatypes. it can be used to improve existing assemblies or in combination with existing automated programs during de novo assembly. the improved version of the musa acuminata genome assembly is accessible and can be downloaded in the new version of the banana genome hub at http://banana-genome.cirad.fr/. tools are available in command line version on github. most of the options are also available on the south green galaxy platform under scaffhunter and scaffremodler toolboxes. availability of supporting data datasets are available through the banana genome hub and the kb library is deposited on the ena read archive. additional file additional file 1: detailed description of tools and processes used to improve the musa acuminata reference sequence and additional figures. abbreviations bacbacterial artificial chromosome gbsgenotyping by sequencing hmwhigh molecular weight lprlarge insert size paired reads ngsnext generation sequencing tetransposable elements wgdwhole genome duplication wgswhole genome shotgun sequencing competing interests the authors declare that they have no competing interests. authors’ contributions gm, fcb, adh: conceived and designed the study and wrote the manuscript. gm: developed the bioinformatic programs and performed the analysis. ak, aa, jd, ah: produced the sequencing data and the genome irys map. gd, mr, ac, ak, jma, ah, fc: contributed to the analysis and edited the manuscript. adh: coordinated the study. all authors read and approved the final manuscript. the authors thank the diversity arrays technology pty ltd for dartseq genotyping, jan vrána and hana Šimková for preparation of hmw dna, and cgiar research program on roots, tubers and bananas for financial support for sequencing data acquisition. we also thank the south green bioinformatics platform for providing us with computational resources. we thank christophe jenny for providing the pahang segregating population from the cirad research station in guadeloupe, french west indies.
we have developed a modular bioinformatics pipeline to improve genome sequence assemblies, which can handle various types of data. bioinformatics tools developed in this work can be used to improve genome sequence assemblies in other species. finally, a consensus gene annotation was projected on the new assembly from two pre-existing annotations. gbs markers were used to anchor scaffolds to pseudo-molecules with a new bioinformatics approach that avoids the tedious step of marker ordering during genetic map construction. the pipeline was used to improve the draft genome sequence of musa acuminata. the release of the musa acuminata reference genome version provides a platform for detailed analysis of banana genome variation, function and evolution. musa acuminatagenome assemblybioinformatics toolpaired-end sequencesgbsgenome mapcirad crp-rtbdart national program of sustainability ilo1204doležel jaroslav issue-copyright-statement© the author 2016 furthermore, a genome map was constructed and used to assemble scaffolds into super scaffolds. long insert size paired-end reads identified scaffold junctions and fusions missed by automated assembly methods. the pipeline comprises several semi-automated tools.
 recent advances in genomics indicate functional significance of a majority of genome sequences and their long range interactions. as a detailed examination of genome organization and function requires very high quality genome sequence, the objective of this study was to improve reference genome assembly of banana.  we have developed a modular bioinformatics pipeline to improve genome sequence assemblies, which can handle various types of data. the pipeline comprises several semi-automated tools. however, unlike classical automated tools that are based on global parameters, the semi-automated tools proposed an expert mode for a user who can decide on suggested improvements through local compromises. the pipeline was used to improve the draft genome sequence of musa acuminata. genotyping by sequencing of a segregating population and paired-end sequencing were used to detect and correct scaffold misassemblies. long insert size paired-end reads identified scaffold junctions and fusions missed by automated assembly methods. gbs markers were used to anchor scaffolds to pseudo-molecules with a new bioinformatics approach that avoids the tedious step of marker ordering during genetic map construction. furthermore, a genome map was constructed and used to assemble scaffolds into super scaffolds. finally, a consensus gene annotation was projected on the new assembly from two pre-existing annotations. this approach reduced the total musa scaffold number from to, with an n that increased from mb to mb.  % of the assembly was anchored to the musa chromosomes compared to the previous 70 %. unknown sites were reduced from to  %.  the release of the musa acuminata reference genome version provides a platform for detailed analysis of banana genome variation, function and evolution. bioinformatics tools developed in this work can be used to improve genome sequence assemblies in other species. electronic supplementary material the online version of this article contains supplementary material, which is available to authorized users. keywords musa acuminatagenome assemblybioinformatics toolpaired-end sequencesgbsgenome mapcirad crp-rtbdart national program of sustainability ilo1204doležel jaroslav issue-copyright-statement© the author 2016 