dataesr/scientific-paragraphs-categorization huggingface dataset:
The results for different sizes of deletions are given in Table  1. Deletions with size between 20 and 250 bp account for a large majority (66. 67%) of all Venter deletions. In this range, JRA clearly outperforms Delly and Lumpy, detecting at least 2. 5 times more true deletions. It is well known that discordant alignments are not very informative in this range, but it is surprising that they can be so detrimental, as witnessed by the drastic loss in sensitivity. This could be attributed to the fact that Delly, for instance, uses its pairedend module and split alignment module in a cascading fashion, such that split alignment is performed only on the regions that is detected by the paired-end module to likely contain SVs. Platypus, in its default mode, detects short deletions (<50 bp) with a very good accuracy, on par with JRA, but with half of its sensitivity.
Evaluation on real sequencing data For evaluation on the CHM1 dataset, we compared the predictions of JRA with Last and the tools used for variant detection: Delly, Lumpy and Platypus, which were provided with the alignments computed by BWA. We omitted the other aligners Splazers and Segemehl because the former failed on larger datasets and the latter required too much running time (see Computational costs). First, comparing the length distribution between the deletions annotated for CHM1  (24)  and for the C. Venter genome (8) (Figure  4A  and B), we can observe that globally, both annotations strongly agree. However, strikingly the shorter deletions in the range 20-35 bp are not reported at all in CHM1. We therefore limited ourselves to predicted deletions which are between 36 and 20 000 bp for all tools.
This observation, along with the discrepancy in the deletion length distributions suggest that the CHM1 annotation might be missing many deletions. To further validate the behavior of our framework on other real genomic data, we used reads sampled from sequencing dataset (SRX1567556) of the Ler-0 strain of the plant A. thaliana, and performed deletion calling with respect to the A. thaliana reference genome. We monitored the number of deletions called as a function of the confidence of the call, and observed similar behavior as in the simulated case (Supplementary Figure  S14 ). For LAST, first pass refers to gapped alignment. For Delly, Lumpy and Platypus, 1st pass refers to BWA alignment, and 2nd pass refers to variant calling including tasks such as sorting of BWA alignments, discordant alignment extraction.
Computational costs We examined running time and memory usage behavior of the four tools using reads from the SRX1567556 A. thaliana dataset and the ERR037900 human whole-genome dataset. From the former, we constructed small datasets by sampling 1 through 5 million reads, in order to understand the growth rate of running time and memory usage. From the latter, we sampled reads to construct datasets of coverages 5, 10, 15 and 20, in order to investigate the performance of the tools on larger-scale input data. All experiments were conducted on a machine with Intel(R) Core(TM) i7-4710MQ CPU and 32GB main memory running Ubuntu 14. 04, and time measurements are for single-threaded operations. Results are shown in Figure  5. Overall running time for JRA consists of time taken by LAST for performing preliminary gapless alignments and time taken for joint alignment.
We can observe that it grows linearly with input size and that it compares well against other tools. It also worth noting that the joint alignment phase imposes only a small overhead. Memory usage of the joint alignment phase depends on the number of reads involved in candidate SV sites, and we can see that it grows linearly with the number of reads. For large-size genomes, we observe that peak memory usage occurs during the preliminary alignment phase. Overall peak memory usage compares well against other tools. Together with the results from performance evaluation, reasonable computational costs further justify the use of JRA as tool for detecting large-sized deletions, even from large datasets.
DISCUSSION We have provided a general framework to incorporate information from a set of error-prone short reads covering an SV site, in order to compute an alignment of this SV region to the reference, while simultaneously reconstructing the SV sequence. We have demonstrated a proof of concept of our joint alignment strategy by applying it to the problem of identifying large deletions and showing its advantage over independent pairwise alignments and over calling tools that rely on such alignments. Given that almost all sequence datasets today are of this nature, it is surprising that current methods rely on pairwise alignments, and that no grouped alignment approach has been proposed before. We hope to have initiated a step in this direction. It is worth noting that in line with our preliminary observation on SV ambiguity (SI-S1), none of the aligners managed to recover more than half of the deletions annotated in the Venter genome, even when considering flawless sequencing.
These are key steps in the current JRA, and to allow complicated evolutionary transformations such as inversions and duplications would require generalizing these steps (See SI-S4 and SI-S3 for further discussions). These are interesting challenges. Another area of improvement is the profile construction step, which currently does not tolerate a high indel error rate. By borrowing ideas from the rich literature in multiple sequence alignment, the scope of application can be expanded to reads with high indel error rates. Another application of our idea of joint alignment is to long read datasets obtained from PacBio or Oxford Nanopore technologies. Long reads allow capturing large variants more accurately, and it also has the added benefit that a single read may span multiple such variants.
Figure5. Growth of running time and memory usage of different methods for down-sampled A. thaliana dataset containing 1 through 5 million reads (top) and human whole genome dataset containing reads corresponding to coverage 5 through 20 (bottom). For JRA, first pass refers to LAST gapless alignments, 2nd pass refers to our joint alignment phase. For LAST, first pass refers to gapped alignment. For Delly, Lumpy and Platypus, 1st pass refers to BWA alignment, and 2nd pass refers to variant calling including tasks such as sorting of BWA alignments, discordant alignment extraction. For Splazers, the first pass alignment is performed by Razers, which gives reads that it cannot fully align to Splazers for the second alignment phase. We did not run Splazers and Segemehl on the whole genome dataset because of insufficient memory for the former and extremely long running time for the latter.
An overlap in the genetic control of T1R and IBD Certain SNP alleles associated with T1R-risk had previously been shown to be susceptibility factors for inflammatory bowel disease (IBD)  [21] [22] [23]. To investigate if there was an enrichment of risk alleles between T1R and IBD, we systematically compared evidence of association with T1R in the Vietnamese discovery set with evidence for association in a recent GWAS meta-analysis for IBD  [24]  (Fig  3A ). Of 232 independent top SNPs that had been associated with IBD by meta-analysis, 208 were available in the T1R-affected and T1R-free GWAS datasets  [24]. For 22/208 SNPs (10. 6%) the IBD risk allele was associated at the 0. 05 level with risk of T1R/leprosy. (Fig 3A  , S6 Table). This observed proportion of shared risk-alleles between T1R leprosy and IBD is significantly non-random (p = 2.
Imputation A total of 38,753 genotyped A/T and C/G SNPs were removed prior to the phasing and imputation. The remaining 495,973 SNPs that passed the quality control filtering in the discovery phase were used to impute additional 11. 5 million variants (SNPs and INDELs) in both T1Raffected and T1R-free family sets with SHAPEIT2  [39]  and IMPUTE2  [40]  software and the 1000 genomes Phase I v3 dataset containing 1092 individuals as the reference panel. Given the exploratory nature of the discovery phase, the threshold of imputation information measure (Info) > 0. 5 was applied to capture most of the common variants (MAF > 5%) with reasonable confidence (S5 Fig)  [41] , MAF > 0. 001 and more than 10 informative families in both T1Raffected and T1R-free sets were used as a post-imputation quality control filtering for the association analyses.
Population-based association analyses were performed using logistic regression under a log-additive model and adjusting by the co-variables gender and age at leprosy diagnosis using PLINK v1. 0. 7. The one-sided test was used with the alternative hypothesis that the T1R-risk alleles were also risk factors in the replication and validation samples. Multivariable analysis were performed with stepwise conditional logistic regression in SAS 9. 3. The haplotype analysis in the Brazilian sample was performed with THESIS v3. 1  [43]. The linkage disequilibrium structure was evaluated with Haploview 4. 1  [44]. To summarize the different steps of the study we used an inverse variance-weighted metaanalysis with a fixed-effect model and an alternative random-effect model proposed by Han and Eskin as implemented in the software METAL  [45]  and METASOFT  [46] , respectively.
For instance, out of the 6,333,954 variants tested for association in our study 319,671 had p < 0. 05 in the T1R-affected subset. Using the observed prior information of the number of variants with p < 0. 05, the hypergeometric test calculates the statistical significance of randomly selecting 22 variants with p < 0. 05 when 208 variants (number of variants from the IBD GWAS meta-analysis present in the T1R dataset) were randomly drawn from a total of ~6. 3 million. Here, the hypergeometric test corresponds to the one-tailed Fisher's exact test. The same analytical approach was applied for the T1R specific variant in IBD, but in this analysis we used the number of variants with p < 0. 05 and p heterogeneity < 0. 05 out of a total of ~6. 3 million variants of the GWAS. Since we tested for sharing of the same risk allele between T1R, and IBD, CD or Ulcerative Colitis (UC) one-tailed p values are reported.
Functional databases and datasets IBD meta-analysis data was freely available at the IBDgenetics website (https://www. ibd genetics. org/) [24, 48]. Briefly, seven CD and eight UC collections with genome-wide data were combined with additional replication samples resulting in a total of 42,950 IBD cases and 53,536 health controls for the IBD meta-analysis  [24, 48]. Variants that surpassed p < 5. 0 x 10 -8 for association with IBD were reported as significant. Functional data for annotated SNPs were extracted from the GTeX (http://www. gtexportal. org/home/) and Haploreg v4 http:// www. broadinstitute. org/mammals/haploreg/haploreg. php databases. [19, 20].
Library preparation and next-generation sequencing For genome assembly, genomic DNA was extracted from a mixed-stage growing population using reagents from the Puregene Core Kit A (QIAGEN, Valencia, CA). Data were generated from a 400-bp library and a 3-kb mate-pair or jumping library, following manufacturer's instructions. The mate-pair library was constructed by the Centre for Genome Research, Liverpool. The 400-bp library was sequenced on the Illumina MiSeq platform (6. 5 million 100-base read pairs and 20. 7 million 300-bp read pairs). The paired-end library was sequenced on an Illumina HiSeq2500 (145. 5 million 100-bp paired reads) (Table  S1 ). Raw data have been submitted to the International Nucleotide Sequence Database Consortium (INSDC) under project accession no. PRJEB15512. RNA was extracted from a mixed-stage growing population of O.
tipulae CEW1, cultured in standard laboratory conditions. Poly(A)-enriched complementary DNA was prepared from the RNA by GATC (Konstanz, Germany) and normalized using reassociation kinetics. RNA sequencing was performed by GATC on the Roche GS FLX Titanium platform. A total of 592,650 reads (average length 369 bp) remained after filtering for quality. Raw transcriptome data are available in INSDC under project no. PRJEB15512). For mapping-by-sequencing, genomic libraries and sequencing data from O. tipulae JU170 (INSDC project accession no. PRJEB19969) and mutant F 2 pools (see project accessions in Table  S3 ) used were generated by BGI. Short insert libraries (,800 bp) were paired-end sequenced on Illumina Hiseq2000, Hiseq2500, or Hiseq4000 with 100-bp reads to obtain 2. 2 GB (403 coverage) of clean data per sample after manufacturer's data filtering (removing adaptor sequences, contamination, and low-quality reads).
Genome assembly All software tools used (including versioning and command line main options) are summarized in Table  S4. Raw reads were trimmed for adaptors using Cutadapt  (Martin 2011)  and low-quality bases, then corrected for sequencing errors based on k-mer content using Quake  (Kelley et al. 2010)  and JELLYFISH  (Marçais and Kingsford 2011). Raw data were checked with FastQC (Andrews 2010) and a preliminary assembly generated with CLC Assembly Cell (CLC bio 2017) (Table  S5 ). The CLC assembly was screened for contaminants using taxon-annotated, GC-coverage (TAGC) plots  (Kumar et al. 2013). Only data deriving from E. coli, the food source, was identified as contaminant and the corresponding reads were removed. The optimal k-mer size for assembly of the cleaned read set was estimated using KmerGenie  (Chikhi and Medvedev 2014).
Gene prediction and orthogroup inference Genes were predicted using a two-pass pipeline  (Koutsovoulos et al. 2014 ) (see Figure  S3 ) based on MAKER2  (Holt and Yandell 2011)  and Augustus  (Stanke et al. 2006) , and using the transcriptome data as evidence. Repeats were identified in the assembly using RepeatModeler  (Smit et al. 2013 (Smit et al. -2015)). MAKER2 was run in an SGE cluster using the SNAP  (Korf 2004 ) ab initio gene finder trained by CEGMA  (Parra et al. 2007 ) output models, the GeneMark-ES ab initio finder, SwissProt proteins, and O. tipulae transcripts. The transcriptome data were filtered so that only reads. 300 bases that had significant basic-localalignment-search-tool (BLAST) similarity to C. elegans protein databases were kept. The MAKER2 predictions were used to train Augustus to generate a custom gene-finding profile for O.
tipulae. Finally, Augustus was used with the gene-finder profile and O. tipulae transcripts to predict the final gene set. Not enough transcript evidence was available to train a model of untranslated regions (UTRs), and therefore no UTRs were annotated. Protein sets for C. elegans (C. elegans Sequencing Consortium 1998), Dictyocaulus viviparus  (Koutsovoulos et al. 2014) , Meloidogyne hapla  (Opperman et al. 2008) , and P. pacificus  (Dieterich et al. 2008) , downloaded from WormBase (http://www. wormbase. org/), were clustered with Orthofinder (Emms and Kelly 2015) using an inflation value of three.
Synteny O. tipulae scaffolds containing. 100 predicted proteincoding genes were selected to perform synteny analyses between O. tipulae and C. elegans. Predicted proteins from these 36 scaffolds were compared to the C. elegans protein data set with BLAST to identify orthologous pairs. For each pair, the chromosome location of the C. elegans ortholog was identified. Hierarchical clustering was performed to group the scaffolds into groups based on the proportions of C. elegans chromosomal attributions (Figure  S9 ).
coli food was just exhausted. Nematodes were directly pooled from different plates and washed in M9 buffer. DNA was then extracted using the Puregene Core Kit A (QIAGEN). Variant analysis, gene mapping, and identification JU170 whole-genome sequencing data were analyzed to identify SNPs compared to the CEW1 reference genome. These variants were then used for genetic mapping of the mutants (listed in Table  S3 ). Reads were mapped with bwa  (Li and Durbin 2009)  to the CEW1 assembly, mappings processed with the GATK tool suite  (McKenna et al. 2010)  version 3. 3-0 and variants called with HaplotypeCaller using default parameters. We followed HaplotypeCaller's authors' recommendations of best practice  (DePristo et al. 2011; Van der Auwera et al. 2013) , realigning reads around indels and performing BQSR by bootstrapping a first call made with HaplotypeCaller.
We analyzed the 300-bp CEW1 MiSeq data used for genome assembly with the same pipeline, after E. coli decontamination, as a control for variant calling. We then hard filtered a list of high-confidence SNPs of JU170 with conservative criteria, retaining homozygous positions covered by at least three reads in each strain, with a sequencing and mapping quality higher than 100 and 55, respectively, and a position noted as reference in CEW1 and variant in JU170. Sequencing data from pooled F 2 mutants were analyzed with the same pipeline, except that variant calling was restricted to a list of JU170 SNPs previously established for faster computing (using the HaplotypeCaller option genotyping_mode GENOTYPE_GIVEN_ALLELES). Output VCF files were used to compute allele frequencies for each SNP on the JU170 list as the ratio of the number of reads with the JU170 allele over the total number of reads.
X chromosome linkage of scaffolds and pyrosequencing Linkage of scaffolds to the X chromosome was determined using directed pyrosequencing of F 1 males from crosses (in both directions) between the CEW1 and JU170 strains. For each scaffold, we selected one polymorphic nucleotide position in the middle of the scaffold in an otherwise conserved context (no other variations in the 300 bp surrounding the variant) to ensure unbiased PCR. For each position, two PCR primers and one sequencing primer were designed (Table  S6 ) using the pyrosequencer's companion design software (PyroMark Q96 ID instrument from Biotage, Uppsala, Sweden). The primers were tested on parental strains. PCR using universal biotinylated primers and single-stranded PCR amplicon purification was performed as previously described  (Duveau and Félix 2010).
Scaffold linkage analysis For each mutant strain (see Table  S3 ) and each scaffold, the mean frequency of alleles in the F 2 mapping population wholegenome sequencing was extracted from the previous pipeline. All data sets included frequencies for the 144 scaffolds (over a total of 191 in nOt. 2. 0) that contained polymorphic positions between the strains CEW1 and JU170. Scaffolds and F 2 mapping populations were clustered using the "heatmap. 2" function of the "gplots" package in R.
Data availability All raw sequencing data supporting the conclusions of this article have been submitted to the INSDC. Sequencing data of genomic DNA and RNA from reference strain CEW1 (used for assembly) are available under project accession no. PRJEB15512. Wholegenome resequencing data of the mapping strain JU170 is under project accession no. PRJEB19969. Project accession numbers corresponding to the sequencing data of F 2 mapping populations of the different strains used in this study are listed in Table  S3. Annotations have also been submitted to the database European Molecular Biology Laboratory-European Bioinformatics Institute under the accession no. PRJEB15512. All scripts used for mapping-by-sequencing and gene identification are available in https://github. com/fabfabBesnard/ Andalusian_Mapping.
Acknowledgments We thank  Sujai Kumar  for a preliminary analysis of the Oscheius tipulae genome. We thank the staff of  Edinburgh Genomics  and the  Centre for Genome Research  for support. We thank  Amhed Vargas Velazquez  for computer help and for pointing out a duplication in the assembly. We thank  Lewis Stevens  for adding O. tipulae nOt. 2. 0 to http:// ensembl. caenorhabditis. org. We thank two anonymous reviewers for their comments and suggestions. We acknowledge WormBase. We gratefully acknowledge support from the  Pôle Scientifique de Modélisation Numérique computing center of École Normale Supérieure de Lyon. G. K. was supported by a  Biotechnological and Biological Sciences Research Council Ph. D. studentship. This work was funded by grants from the  Agence Nationale de la Recherche  ( ANR12-BSV2-0004-01  and  ANR10-LABX-54   MEMOLIFE ).
A standardized and scalable approach The EMP solicited the global scientific community for environmen tal samples and associated metadata spanning diverse environments and capturing spatial, temporal, and/or physicochemical covariation. The first 27,751 samples from 97 independent studies (Supplementary Table  1 ) represent diverse environment types (Fig. 1a ), geographies (Fig. 1b ), and chemistries (Extended Data Fig. 1 ). The EMP encom passes studies of bacterial, archaeal, and eukaryotic microbial diversity. The analysis here focuses exclusively on the bacterial and archaeal components of the overall database (for concision, use of 'microbial' will hereafter refer to bacteria and archaea only). Associated meta data included environment type, location information, host taxonomy (if relevant), and physico chemical measurements (Supplementary Table  2 ).
Physicochemical measurements were made in situ at the time of sampling. Investigators were encouraged to measure temperature and pH at minimum. Salinity, oxygen, and inorganic nutrients were measured when possible, and investigators collected additional meta data pertinent to their particular investigations. Metadata were required to conform to the Genomic Standards Consortium's MIxS and Environment Ontology (ENVO) standards  10, 11. We also used a lightweight application ontology built on top of ENVO: the EMP Ontology (EMPO) of microbial environments. EMPO was tailored to capture two major environmental axes along which micro bial betadiversity has been shown to orient: host association and salinity 1,2. We indexed the classes in this application ontology (Fig. 1a ) as levels of a structured categorical variable to classify EMP samples as hostassociated or freeliving (level 1).
We generated sequence reads of 90-151 base pairs (bp) (Extended Data Fig. 2a , Supplementary Table  1 ), totaling 2. 2 billion sequences, an average of 80,000 sequences per sample. Sequence analysis and taxonomic profiling were done initially using the common approach of assigning sequences to operational taxonomic units (OTUs) clustered by sequence similarity to existing rRNA data bases  14, 20. While this approach was useful for certain analyses, for many sample types, especially plantassociated and freeliving communities, onethird of reads or more could not be mapped to existing rRNA databases (Extended Data Fig. 2b ). We therefore used a referencefree method, Deblur  21  , to remove suspected error sequences and provide singlenucleotide resolution 'subOTUs' , also known as 'amplicon sequence variants'  22  , here called 'tag sequences' or simply 'sequences'.
Microbial ecology without OTU clustering While earlier largescale 16S rRNA amplicon studies adopted OTU clustering approaches in part out of concern that erroneous reads would dominate diversity assessments  23  , patterns of prevalence (presenceabsence) in our results suggest that Deblur error removal produced ecologically informative sequences without clustering. After rarefying to 5,000 sequences per sample, a total of 307,572 unique sequences were contained in the 96 studies and 23,828 samples of the 'QCfiltered' Deblur 90bp observation table. Among studies, more than half (57%) of all obtained sequences were observed in two or more studies, but only 5% were observed in more than ten studies; the most prevalent sequence was found in 88 of 96 studies (Extended Data Fig. 4a ). Among samples, although most sequences (86%) were observed in two or more samples, only 7% were observed in more than 100 samples (Extended Data Fig.