Robert F Rmch

Msc Julia ; Lee

King'

MD Bronwen E Shaw
email: bshaw@doctors.org.uk

MD J F Apperley

MD N H Russell

MD C Craddock

MD E Liakopoulou

MD M N Potter

MD R Wynn

MD B Gibson

Msc M Pearce

Msc K Kirkland

MD J A Madrigal

MD G Cook

MD J L Byrne

R Pearce

Unrelated donor peripheral blood stem cell transplants incorporating pretransplant in-vivo Alemtuzumab are not associated with any increased risk of significant acute or chronic GVHD

Keywords: unrelated donor, T cell depletion, myeloablative, PBSC, BM

Unrelated donor peripheral blood stem cell transplants incorporating pre-transplant

Introduction

The use of unrelated donors for allogeneic transplantation has extended the availability of this potentially life-saving therapy to those who lack a sibling donor. In view of the increased genetic disparity between recipient and donor the incidence of some complications, such as graft versus host disease (GvHD), are increased. Recent advances in technology, including those related to HLA-typing, have resulted in the outcome of UD transplants reaching that of transplants using sibling donors (Goldstone andRowe 2009, Lee, et al 2009), in part due to a reduction in the incidence of severe GvHD. One highly successful strategy to reduce GvHD is the use of T cell depletion (TCD) -either in vivo or ex vivo [START_REF] Kottaridis | In vivo CAMPATH-1H prevents graftversus-host disease following nonmyeloablative stem cell transplantation[END_REF][START_REF] Marks | T cell-depleted unrelated donor bone marrow transplantation for acute myeloid leukemia[END_REF]. A commonly used method of in-vivo TCD is Alemtuzumab, given for several days pretransplant. Alemtuzumab is an antibody against CD52 and as such depletes not only T-cells, but also other cell subset such as B cells and dendritic cells [START_REF] Hale | Improving the outcome of bone marrow transplantation by using CD52 monoclonal antibodies to prevent graft-versus-host disease and graft rejection[END_REF]. Alemtuzumab has been shown to persist in the patient's circulation for up to 30 days post infusion, allowing for depletion of CD52-positive cells from the patient as well as from the donor cells on infusion [START_REF] Rebello | Pharmacokinetics of CAMPATH-1H in BMT patients[END_REF]. Historically all donors provided stem cells via a bone marrow (BM) harvest. Since the 1990s we have had the ability to collect GCSF -mobilised peripheral blood stem cells (PBSC) and since 2000 their use in UD has increased considerably [START_REF] Ljungman | Allogeneic and autologous transplantation for haematological diseases, solid tumours and immune disorders: current practice in Europe 2009[END_REF]. In view of the low incidence of significant short-term side effects and the ease of donation, approximately 70-80% of unrelated donors opt to donate via this method in preference to BM harvesting [START_REF] Miller | Recovery and safety profiles of marrow and PBSC donors: experience of the National Marrow Donor Program[END_REF][START_REF] Pulsipher | Adverse events among 2408 unrelated donors of peripheral blood stem cells: results of a prospective trial from the National Marrow Donor Program[END_REF]. It was clear from early experience that the number of cells (CD34+ and CD3+ cells) harvested via PBSC was greater than that by BM harvesting [START_REF] Pavletic | Lymphocyte reconstitution after allogeneic blood stem cell transplantation for hematologic malignancies[END_REF]. While this was associated with a significantly shorter engraftment times, there were concerns that the incidence of GvHD would be increased with the use of PBSC. Indeed, a number of studies, predominantly in sibling transplants, showed an increased incidence of both acute and/or chronic GvHD in PBSC recipients [START_REF] Schmitz | Filgrastim-mobilized peripheral blood progenitor cells versus bone marrow transplantation for treating leukemia: 3-year results from the EBMT randomized trial[END_REF][START_REF] Schmitz | Long-term outcome of patients given transplants of mobilized blood or bone marrow: A report from the International Bone Marrow Transplant Registry and the European Group for Blood and Marrow Transplantation[END_REF], Storek, et al 1997) without, however, an impact on patient outcome (Friedrichs, et al , Gallardo, et al 2009[START_REF] Ringden | Transplantation of peripheral blood stem cells as compared with bone marrow from HLAidentical siblings in adult patients with acute myeloid leukemia and acute lymphoblastic leukemia[END_REF][START_REF] Schmitz | Filgrastim-mobilized peripheral blood progenitor cells versus bone marrow transplantation for treating leukemia: 3-year results from the EBMT randomized trial[END_REF]. More recent studies have considered UD transplants [START_REF] Blau | Feasibility and safety of peripheral blood stem cell transplantation from unrelated donors: results of a single-center study[END_REF], Eapen, et al 2007, Garderet, et al 2003[START_REF] Remberger | Increased risk of extensive chronic graft-versus-host disease after allogeneic peripheral blood stem cell transplantation using unrelated donors[END_REF][START_REF] Remberger | No difference in graft-versus-host disease, relapse, and survival comparing peripheral stem cells to bone marrow using unrelated donors[END_REF], again confirming an increase in GvHD with PBSC, but not consistently showing any difference in outcome. Most of these studies, however include recipients of both HLA matched and mismatched grafts and all report series where transplants are predominantly T cell replete. A few studies have suggested a differential effect on survival due to stem cell source, both in sibling (Eapen, et al 2004[START_REF] Group | Allogeneic peripheral blood stem-cell compared with bone marrow transplantation in the management of hematologic malignancies: an individual patient data meta-analysis of nine randomized trials[END_REF] and UD transplants (Eapen, et al 2007, Garderet, et al 2003). The large meta-analysis including 1111 adult recipients of HLA-identical sibling grafts [START_REF] Group | Allogeneic peripheral blood stem-cell compared with bone marrow transplantation in the management of hematologic malignancies: an individual patient data meta-analysis of nine randomized trials[END_REF], showed that PBSC was associated with increased incidences of both acute and chronic GvHD. However, PBSC was also associated with a decrease in relapse in both early and late stage disease, with a significantly improved survival in recipients of PBSC with late stage disease. Conversely, in a later study (Eapen, et al 2004) in children receiving HLA-identical sibling allografts, PBSC was associated with a higher incidence of chronic GvHD, as well as a worse overall survival and a higher transplant related mortality (TRM). Furthermore in a CIBMTR study of UD transplant outcome, Eapen et al (Eapen, et al 2007) showed a significantly higher incidence of both grade II-IV acute GvHD and of chronic GVHD. Although there was no detrimental effect on overall survival in the whole group, the subgroup of patients with chronic myeloid leukaemia (CML) in first chronic phase (CP1) receiving PBSC did have a significantly worse overall survival, with a trend to an increase in TRM. In the study by Gardaret et al (Garderet, et al 2003), there was a trend to an increase in acute (but not chronic) GvHD in patients with ALL receiving PBSC, which was associated with an inferior survival compared to BM (without an impact on relapse). We hypothesised that any negative impact caused by increased GvHD using PBSC compared to BM might be abrogated in recipients of UD transplants incorporating in-vivo Alemtuzumab conditioning for T cell depletion. Here we present data on the outcome in 306 UK transplant recipients with leukaemia studying the impact of stem cell source on outcome.

Materials and methods

Data collection

Patients fulfilling the requirements for this study were identified from the database held by the British Society for Blood and Marrow Transplantation (BSBMT) registry in London, UK. Data are reported at fixed time points in the form of the MED A and B forms derived from the European Group for Blood and Marrow Transplant (EBMT)). Data collected by the Anthony Nolan Trust from a previous study were added to the BSBMT database. 306 patients from 19 centres were identified. Supplementary data requests were sent to all centres. Criteria for inclusion were: 1. a transplant at a UK centre between 2000 and 2007. 2. myeloablative conditioning 3. in-vivo pre-transplant T cell depletion using alemtuzumab 4. HLA matched (10/10) unrelated donor 5. standard risk leukaemia (acute myeloid or lymphoblastic leukaemia in first or second complete remission (CR), CML in CP1).

Patients, definitions and inclusion criteria

HLA typing was performed in the local Histocompatibility and Immunogenetics laboratories by molecular methods to achieve high resolution types.

Alemtuzumab schedule and dose were given according to institutional protocols. Both the overall dosing and the days of infusion differed between institutions. However the dose of Alemtuzumab within an institution did not differ depending on whether the stem cell product was BM or PBSC. CMV was treated pre-emptively. Donor leukocyte infusions (DLI) were not used routinely, but given in some cases according to institutional protocols for mixed chimerism, relapse or viral reactivations.

Primary Graft Failure was defined as a failure to achieve a neutrophil count of 0.5 x 10 9 /l by day 28 and was evaluable only in those surviving at least 28 days after the transplant. Relapses were defined as haematological, cytogenetic or molecular. Acute and chronic GvHD were defined using internationally accepted criteria.

The protocols used were approved by the individual institutional review boards of all of the contributing hospitals. Informed consent for the transfer of data to, and the analysis by, the BSBMT was obtained from all patients.

Statistical analyses

Relationships between categorical variables were analysed by Fisher's exact test. Continuous variables such as age were tested using the k-sample equality of medians test. Overall survival was calculated by Kaplan-Meier analysis. Univariate analysis of OS and PFS was performed using the log-rank test, and multivariate analysis by Cox regression. Cumulative incidence of non-relapse mortality was calculated by competing risks regression, with death from relapse as the competing risk. Relapse rate was calculated as a cumulative incidence by competing risks regression, with transplant related death and chronic GvHD as the competing risks. Cumulative incidence of cGvHD was calculated by competing risks regression with death and relapse as competing risks. All multivariate models were designed to include any factor which had a p value of <0.2 in univariate analysis plus the source of stem cells. In view of the relatively small dataset and number of significant variables, stepwise inclusion / exclusion were not applied.

Statistics tests were performed using Stata (StataCorp, Texas 77845 URL: http://www.stata.com), with competing risks calculated using "stcompet" (Enzo Coviello, Italy; May Boggess, StataCorp); and using R (R Development Core Team (2006). R Foundation for Statistical Computing, Vienna, Austria. ISBN 3-900051-07-0, URL http://www.R-project.org), with competing risks calculated using the package "cmprsk" (RJ Gray). Patient and transplant characteristics 306 patients were included in the study. Patient, disease and transplant characteristics are presented in Table 1. The majority of conditioning regimens included TBI (92%), usually in combination with cyclophosphamide (84%), and additional GvHD prophylaxis with cyclosporine alone or cyclosporine and short-course methotrexate was used. In all patients T cell depletion was performed using pre-transplant in vivo Alemtuzumab. The two groups were similar, except for an increased use of BM in patients with CML (p=0.007) and the disease stage differed significantly between the groups (p=0.001). 24 patients received donor lymphocyte infusions (DLI) posttransplant, 17 for disease relapse, 2 for viral reactivations and 5 for mixed chimerism. Engraftment 1% (4/306) of patients had primary graft failure. Of these, 2 received BM and 2 received PBSC (two had a second allograft and are alive, while two succumbed from infection or disease progression, one each in the BM and PBSC groups). A further 4 patients who received BM had secondary graft failure (two of whom remain alive following further cellular therapy) while none of the PBSC recipients developed this complication. There was no statistically significant difference between stem cell sources in the incidence of graft failure overall (p=0.306). The median time to neutrophil engraftment was 18 days (range: 9-61 days). This was significantly faster in recipients of PBSC (15 days) than BM (20 days, p<0.0001). Likewise, platelet engraftment was significantly faster in recipients of PBSC (17 days) than BM (28 days, p<0.0001), with a median of 24 days (5-469).

Graft versus Host Disease (GvHD)

All patients who achieved neutrophil engraftment were considered eligible for analysis. The incidence of acute GvHD (aGvHD) of any grade, grade II-IV and grade III-IV were 56%, 23% and 4% respectively. There was a significantly higher risk of developing aGvHD in the recipients of PBSC (65%) compared to BM recipients (49%, p=0.012), but this risk was only apparent in grade 1 disease. Most importantly, there were no significant differences in severe aGvHD between the groups (grade III/IV: BM 4% PBSC 5%, p=0.554) (Table 2). No other factor had a significant impact on the incidence or severity of aGvHD in either univariate or multivariate analysis.

Overall chronic GvHD (cGvHD) developed in 119/235 patients (48% at 5 years). There was no significant difference in the incidence of cGvHD dependant on stem cell source (5 years: 47% following BMT and 49% following PBSCT; p=0.48) (figure 1). Strikingly, extensive cGvHD occurred in only 15% of BM recipients and 13% of PBSC recipients (p=0.831). Factors associated with cGvHD were CMV seropositivity in the patient (Hazard ratio (HR) =1.56; Confidence interval (CI) 1.07 -2.27, p=0.029) and disease in CR1/CP1 (HR 0.66; CI 0.43 -1.00, p=0.054). In multivariate analysis, including these two factors and source of stem cells and age, only disease status at transplant showed a trend towards significance (HR=0.68; CI 0.44 -1.05, p=0.081).

Non-relapse mortality (NRM)

NRM was 9% and 20% at day 100 and 1 year respectively (figure 2). There were no significant differences in NRM between recipients of BM or PBSC (day 100: 9% and 11%, 1 year: 20% and 20%, and 5 year: 23% and 25% respectively, p=0.77). There was no significant impact on NRM 7 due to the presence or absence of aGvHD. Factors associated with NRM in univariate analysis were: patient CMV (HR 2.15; CI 1.34 -3.47, p=0.001) and donor (HR 1.54; CI 0.93 -2.53, p=0.090) CMV seropositivity. In multivariate analysis (including these factors, patient age and stem cell source) the only factor to remain significantly associated with a higher mortality was CMV seropositivity in the patient (HR=1.98; CI 1. 16 -3.38, p=0.012). The causes of NRM were: GvHD (22%), infection (62%) and other (16%), which did not differ significantly depending on stem cell source (p=0.358).

Relapse

The cumulative incidences of relapse were 36% and 41% at 2 and 5 years respectively (figure 3). Relapse was predominantly haematological in patients with acute leukaemia (86%), with 6% defined as cytogenetic and 8% defined as molecular. In contrast, only 41% of relapses in CML were haematological with 26% defined as cytogenetic and 33% defined as molecular. There was a trend to a difference in relapse rate according to whether the patient received BM or PBSC (5 years: 45% and 33% respectively, p=0.096). The only factor to impact significantly on disease relapse was disease type, with a relapse risk of 58% at 5 years in CML patients compared to 33% in AL patients (p<0.001). There was a trend to an increase in relapse in CMV seronegative patients (HR=0.69; CI 0.46-1.04, p=0.075). There were no differences in relapse risk between patients with AL transplanted in CR1 or CR2 (both 33%). In multivariate analysis, including disease, patient and donor CMV and stem cell source, the only significant factor was the underlying disease (HR=2.03; CI 1.38-2.98, p<0.001). There was no significant impact on relapse due to the presence or absence of aGvHD.

Survival

The median follow up is 48 months (3-110 months). This was significantly longer in recipients of BM than PBSC (59 vs. 36 months, p=0.001) (figure 4). The overall survival in the whole group was 53% at 8 years, with progression free survival of 36% at the same time point. There were no significant differences in survival between patient who received BM or PBSC (8 years: 54% and 52% respectively, p=0.571). Likewise progression free survival was similar (8 years: 32% and 42% respectively, p=0.225). There was no significant impact on OS due to the presence or absence of aGvHD. Two factors were significantly associated with an improved survival in univariate analysis: CML compared to AL (p=0.005) and CMV seronegativity in patients (p=0.025). Both of these factors remained significant in multivariate analysis (HR 0.57; CI 0.37-0.86, p=0.010 and HR 0.70; CI 0.49-0.99, p=0.046 respectively), while PBSC vs. BM was not significant (HR 1.03;p=0.826).

Discussion

This study is the first to consider the impact of stem cell source in a cohort of in-vivo T cell depleted HLA matched UD allograft recipients. In contrast to previous studies of T-replete UD transplants, we show no increase in the risk of grade II-IV acute GVHD or of chronic GVHD in recipients of PBSC. The only significant difference between BM and PBSC recipients in this study was a higher incidence of grade I aGvHD in the PBSC group.

In 617 adult recipients of T cell replete UD transplants, Eapen et al (Eapen, et al 2007) reported an increase in the incidence of overall, but not grade III/IV, aGvHD with PBSC. In addition, F o r P e e r R e v i e w unlike our report, the incidence of cGvHD was also found to be significantly higher in the patients receiving PBSC. Two smaller studies by Gardaret et al (Garderet, et al 2003) (n=213) and Blau et al [START_REF] Blau | Feasibility and safety of peripheral blood stem cell transplantation from unrelated donors: results of a single-center study[END_REF] (n=74) reported significant increases in aGvHD (only in HLA mismatched patients in the latter study), but no increase in cGvHD. In two studies by Rembeger et al, reporting on both short [START_REF] Remberger | No difference in graft-versus-host disease, relapse, and survival comparing peripheral stem cells to bone marrow using unrelated donors[END_REF] and long term [START_REF] Remberger | Increased risk of extensive chronic graft-versus-host disease after allogeneic peripheral blood stem cell transplantation using unrelated donors[END_REF] follow-up in 214 adult transplant recipients, there was no increase in aGvHD while the incidence of chronic extensive GvHD was significantly higher when PBSC was used compared to BM.

There are a number of clear differences between these studies, the most important being the use or not of T-cell depletion, in particular, Alemtuzumab. A possible additional factor explaining the lack of impact of stem cell source on cGvHD is the differing length of follow-up in each study. The median follow-up in the Remberger study is over 4 years (which is the same as in our current study), while that in the other studies ranges between 1.5 and 3 years and, in some, differs between the PBSC and BM groups. The pattern of GvHD may also be of interest. All of the studies to date report outcomes using the 'old' classification of GvHD i.e. using a simply time cut off between acute and chronic. It may be that using the reclassification of GvHD in future studies, would provide further insights (Filipovich 2008). In our study we noticed that the time to onset of cGVHD was shorter in the PBSC than the BM group (data not shown), suggesting that a proportion of the 'chronic' GvHD would now be reclassified as late onset acute or persistent rather than de novo chronic GvHD. It is well recognised that chronic GvHD is more likely in patients who have suffered acute GvHD (especially of higher grades) which may explain the increase in cGvHD in the Eapen study (where the incidence of acute GvHD was higher than in our study), while this may be less discernable in the Blau and Gardaret study where the overall and grade II/IV aGvHD incidence was less.

T cell depletion using Alemtuzumab was universal in our study. This method of GvHD prophylaxis was either absent or used in a very small minority of patients in the other studies except for the Remberger study where one third of patients received ATG. Direct comparisons between ATG and Alemtuzumab cannot be made given that Alemtuzumab depletes all CD52positive cells and not just T cells. The use of ATG in only a proportion of patients and the likelihood of 'hidden' HLA mismatches in that study, especially for HLA-C (low resolution HLA-A and -B and high resolution -DRB1 matching criteria were used) may account for the high incidence of chronic GvHD that they observed. All of the recipients included in our study had high resolution typing performed. The significant reduction in overall incidence of GvHD following TCD strategies (in particular grade III-IV acute and extensive chronic) compared to that following T cell replete transplantation has been reported in both randomised [START_REF] Bacigalupo | Thymoglobulin prevents chronic graft-versus-host disease, chronic lung dysfunction, and late transplant-related mortality: long-term follow-up of a randomized trial in patients undergoing unrelated donor transplantation[END_REF][START_REF] Bacigalupo | Antithymocyte globulin for graft-versus-host disease prophylaxis in transplants from unrelated donors: 2 randomized studies from Gruppo Italiano Trapianti Midollo Osseo (GITMO)[END_REF], Finke, et al 2009, Wagner, et al 2005) as well as non-randomised studies (Chakraverty, et al 2002[START_REF] Malladi | Alemtuzumab markedly reduces chronic GVHD without affecting overall survival in reduced-intensity conditioning sibling allo-SCT for adults with AML[END_REF][START_REF] Kennedy-Nasser | Comparable outcome of alternative donor and matched sibling donor hematopoietic stem cell transplant for children with acute lymphoblastic leukemia in first or second remission using alemtuzumab in a myeloablative conditioning regimen[END_REF][START_REF] Kottaridis | In vivo CAMPATH-1H prevents graftversus-host disease following nonmyeloablative stem cell transplantation[END_REF][START_REF] Malladi | Alemtuzumab markedly reduces chronic GVHD without affecting overall survival in reduced-intensity conditioning sibling allo-SCT for adults with AML[END_REF][START_REF] Mohty | Graft-versus-host disease following allogeneic transplantation from HLA-identical sibling with antithymocyte globulin-based reduced-intensity preparative regimen[END_REF]. Conversely, both infection [START_REF] Bacigalupo | Antithymocyte globulin for graft-versus-host disease prophylaxis in transplants from unrelated donors: 2 randomized studies from Gruppo Italiano Trapianti Midollo Osseo (GITMO)[END_REF][START_REF] Malladi | Alemtuzumab markedly reduces chronic GVHD without affecting overall survival in reduced-intensity conditioning sibling allo-SCT for adults with AML[END_REF] and relapse (Wagner, et al 2005) may be increased in the TCD setting, due to delayed immune reconstitution and a reduced graft versus leukaemia effect. In this study CMV seropositivity in the recipient was associated with an inferior survival. It is well recognised that CMV reactivation is common following TCD (Finke, et al 2009, Wagner, et al 2005). Although the direct impact of CMV reactivation on mortality is minimal due to pre-emptive strategies, it is known that CMV seropositivity is associated with a predisposition to other infections and may impact on survival in other less well defined ways (Craddock, et al 2001[START_REF] Nichols | High risk of death due to bacterial and fungal infection among cytomegalovirus (CMV)-seronegative recipients of stem cell transplants from seropositive donors: evidence for indirect effects of primary CMV infection[END_REF]. Although the relapse risk in the AL patients in this study was acceptable, the relapse risk in CML was disappointingly high (although not significantly different based on stem cell source). This may be due to the inclusion of TCD in the conditioning as previously reported (Wagner, et al 2005).

Neither our or other UD studies have reported any significant differences due to stem cell source in TRM, relapse or OS in the overall group studied. In the Eapen study there was a trend to a higher TRM, with a significantly lower OS, in recipients with CML in CP1, although the reason for this is not clear. In the Gardaret study patients with ALL receiving PBSC had a significantly lower DFS and OS. Besides the obvious differences with regards to TCD, other differences between the studies are noted, such as the inclusion of both paediatric and adult patients, the degree of HLA matching, the diseases and stage of disease included, medical GvHD prophylaxis and the year of transplantation. None of the studies in UDs are randomised and the results of the BMT CTN 0201 randomised trial (Peripheral Blood versus Bone Marrow Grafts from Unrelated Donors) (http://www.cibmtr.org/Studies/ClinicalTrials/BMT_CTN/Protocols/0201/index.html) are keenly awaited.

It is recognised that the current study has a number of limitations. In common with most retrospective studies there is some heterogeneity in the study group. Although TCD using Alemtuzumab was universal, the schedule and dosage varied between centres. Importantly, however, the dose of Alemtuzumab within an institution did not differ depending on whether the stem cell product was BM or PBSC. In addition the use of cyclosporin (with or without methotrexate) differed by centre. DLI were given for various reasons, but none as part of a planned immunotherapy program. Any impact of DLI on outcomes should be negligible, as most were given for relapse, and this was taken into account as a competing risk for chronic GVHD. Although the number of CD34+ cells is known, the data on CD3+ and other cell subset content of the graft were not available.

In conclusion, in this study we confirm the low incidence of severe GvHD following TCD myeloablative transplantation from UD for patients with standard risk leukaemia, with excellent long-term survival rates. The use of either BM or PBSC was associated with equally good outcomes. The only significant difference between the groups was an increase in grade 1 acute GvHD in the PBSC recipients. The risk of extensive cGVHD was low and did not differ with stem cell source. These data suggest that either stem cell source can be used with equivalent results and that there is no deleterious effect from the use of PBSC in Alemtuzumab containing unrelated donor transplants.

  . (1997) Allogeneic peripheral blood stem cell transplantation may be associated with a high risk of chronic graft-versus-host disease. Blood, 90, 4705-4709. Wagner, J.E., Thompson, J.S., Carter, S.L. & Kernan, N.A. (2005) Effect of graft-versus-host disease prophylaxis on 3-year disease-free survival in recipients of unrelated donor bone marrow (T-cell Depletion Trial): a multi-centre, randomised phase II-III trial. Lancet,

Figure 1 :

 1 Figure legends:

Figure 2 :

 2 Figure 2: Non-relapse mortality using either PBSC or BM.

Figure 3 :Figure 4 :

 34 Figure 3: The incidence of relapse using either PBSC or BM

Table 2 :

 2 Incidence of acute GvHD comparing PBSC to BM. Refers to the lack of statistical significance between those with grade III/IV acute GvHD in the BM and PBSC group (i.e grade III/IV vs 0-II).

	Acute GvHD	BM	PBSC	P value
	Present	87 (49%)	77 (65%)	0.012
	Grade I	46 (26%)	48 (40%)
	Grade II	33 (19%)	23 (19%)
	Grade III	5 (3%)	4 (3%)	0.554*
	Grade IV	1 (1%)	2 (2%)
	Unknown grade	2	0
	*

British Journal of Haematology

Acknowledgements

We would like to thank all of the UK transplant teams who have contributed patients to this study.

Contributors

Authorship contributions

JLB, NHR and BES contributed to the design and analysis. RP performed the statistical analysis. All authors contributed to the data collection, writing and review of the manuscript.

Conflict of Interest

There are no conflicts of interest to disclose for any of the authors. Chakraverty, R., Peggs, K., Chopra, R., Milligan, D.W., Kottaridis, P.D., Verfuerth, S., Geary, J., Thuraisundaram, D., Branson, K., Chakrabarti, S., Mahendra, P., Craddock, C., Parker, A., Hunter, A., Hale, G., Waldmann, H., Williams, C.D., Yong, K., Linch, D.C., Goldstone, A.H. & Mackinnon, S. (2002) Limiting transplantation-related mortality following unrelated donor stem cell transplantation by using a nonmyeloablative conditioning regimen. Blood, 99, 1071-1078. Craddock, C., Nagra, S., Peniket, A., Brookes, C., Buckley, L., Nikolousis, E., Duncan, N., Tauro, S., Yin, J., Liakopoulou, E., Kottaridis, P., Snowden, J., Milligan, D., Cook, G 11,[331][332][333][334][335][336][337][338]D.,de la Camara,R.,Nieto,J.B.,Espigado,I.,Iriondo,A.,A.,Vallejo,C.,Martin,C.,Caballero,D.,Brunet,S.,Serrano,D., Solano, C., Ribera, J.M., de la Rubia, J. & Carreras, E. ( 2009) Is mobilized peripheral blood comparable with bone marrow as a source of hematopoietic stem cells for allogeneic transplantation from HLAidentical sibling donors? A case-control study. Haematologica, 94, 1282-1288. Garderet, L., Labopin, M., Gorin, N.C., Polge, E., Fouillard, L., Ehninger, G.E., Ringden, O., M -male, F -female, TBI -total body irradiation, Cy -cyclophosphamide, Bu -busulphan, AML -acute myeloid leukaemia, ALL -acute lymphoblastic leukaemia, CML -chronic myeloid leukaemia, AL -acute leukaemia, CR -complete remission, CP -chronic phase


Eric Rivals

Leena Salmela

Jorma Tarhio

EXACT SEARCH ALGORITHMS FOR BIOLOGICAL SEQUENCES

With the development of sequencing techniques, it has become easy to obtain the sequence, i.e. the linear arrangement of residues (nucleotides or amino-acids), of DNA, RNA, or protein molecules. However, determining the function of a molecule remains difficult and is often bound to find a sequence similarity to another molecule whose role in the cell is at least partially known. Then, the biologist can predict that both molecules share the same function and try to check this experimentally. Functional annotations are transferred from one sequence to the other provided that their similarity is high enough. This procedure is also applied to molecules subparts, whose sequences are shorter: protein domains, DNA/RNA motifs, etc.

Depending on the sequence lengths and expected level of evolutionary relatedness, the sequence similarity can be found using alignment or pattern matching procedures. A quest in bioinformatics has been to design more sensitive sequence similarity searching methods to push further the limit or gray

zone at which evolutionary sequence similarity cannot be departed from random sequence similarity [START_REF] Altschul | Gapped blast and psi-blast: a new generation of protein database search programs[END_REF][START_REF] Gusfield | Algorithms on Strings, Trees and Sequences[END_REF]. These methods (e.g. profile Hidden Markov Models) have provided, at the expense of computing time, important improvements in functional annotations. However, it has soon become clear that in other frameworks only high level similarity was sought, and speed rather than sensitivity was the major issue. Hence, researchers have designed a continuum of methods that can be classified according to level of allowed dissimilarity:

1. full sensitivity alignment (Smith & Waterman algorithm [START_REF] Smith | Identification of common molecular subsequences[END_REF]), 2. fast similarity search programs (e.g. BLAST [START_REF] Altschul | Gapped blast and psi-blast: a new generation of protein database search programs[END_REF]), 3. approximate pattern matching (e.g. Bowtie [START_REF] Langmead | Ultrafast and memory-efficient alignment of short dna sequences to the human genome[END_REF]), 4. near-exact and exact pattern matching (e.g. mpscan [START_REF] Rivals | Mpscan: fast localisation of multiple reads in genomes[END_REF]).

For some of everyday sequence manipulation tasks, the user needs exact pattern matching programs (as available in large bioinformatic program suites like EMBOSS): to find back from which chromosome or where in a genome does a given sequence comes from; to find short nucleotidic motifs, like restriction or cleavage sites, in long DNA sequences; to verify whether a distinguishing sequence motif really departs negative from positive instances (longer sequences). The latter happens when designing oligonucleotides for gene expression arrays or multiple primers for multiplex polymerase chain reactions [START_REF] Rahmann | Fast large scale oligonucleotide selection using the longest common factor approach[END_REF]. Even for exploring protein sequences, a server has been launched that offers an exact search for short polymers in all sequences of protein data-banks [START_REF] Bleasby | OWL -A non-redundant, composite protein sequence database[END_REF]. In such frameworks, the need is single pattern search or multiple pattern search for a few hundreds patterns, which can easily be solved by repetitively applying a single pattern matching program. Algorithmic solutions for these tasks will be described in Section 1.2. However, pattern matching algorithms fail to become popular among biologists for several reasons:

• most of them lack implementations capable of handling biological sequence formats (which then requires to change the format),

• they lack a graphical interface or were not integrated in popular graphical sequence exploration package like the GCG (Genetics Computer Group) package,

• as BLAST [START_REF] Altschul | Gapped blast and psi-blast: a new generation of protein database search programs[END_REF] was used for similarity searching on a daily basis, it has become the all-purpose tool for most sequence processing tasks, even when more adapted solutions were available [START_REF] Claverie | Bioinformatics for dummies[END_REF].

Since 2005, biology experiences the revolution of high-throughput sequencing (HTS) due to the renewal of sequencing techniques (new technologies are often termed Next Generation Sequencing) [START_REF] Margulies | Genome sequencing in microfabricated high-density picolitre reactors[END_REF]. Due to the invention of parallel sequencing of multiple molecules on a single machine, the sequencing output per run has grown by orders of magnitude compared to the traditional Sanger technique, and is expected to increase further [START_REF] Editor | Prepare for the deluge[END_REF]. This change has not only technological consequences. Experiments previously done by hybridisation are now preferentially performed by sequencing [START_REF] Blow | Transcriptomics: The digital generation[END_REF], since these techniques offer a much deeper sampling and allows to cover the whole genome. Hence, HTS is now exploited to address surprisingly diverse biological questions: genome sequencing or resequencing [START_REF] Margulies | Genome sequencing in microfabricated high-density picolitre reactors[END_REF][START_REF] Aury | High quality draft sequences for prokaryotic genomes using a mix of new sequencing technologies[END_REF], transcriptomics [START_REF] Sultan | A Global View of Gene Activity and Alternative Splicing by Deep Sequencing of the Human Transcriptome[END_REF][START_REF] Philippe | Using reads to annotate the genome: influence of length, background distribution, and sequence errors on prediction capacity[END_REF], genomic variation identification or genome breakpoint mapping [START_REF] Chen | Mapping translocation breakpoints by next-generation sequencing[END_REF], metagenomics [START_REF] Huson | Methods for comparative metagenomics[END_REF], and epigenomics [START_REF] Boyle | High-Resolution Mapping and Characterization of Open Chromatin across the Genome[END_REF]. To grasp how drastic the shift is, consider that one epigenomic census assays published in 2007 produced in one experiment an already amazing 1.5 millions short read sequences of 27 base pairs (bp)1 each [START_REF] Johnson | Genome-Wide Mapping of in Vivo Protein-DNA Interactions[END_REF], while another published only one year later delivered with the same technology 15 millions 20 bp reads [START_REF] Boyle | High-Resolution Mapping and Characterization of Open Chromatin across the Genome[END_REF].

In all such applications, the first bioinformatic task is to map the short reads on a reference genome sequence or on a large collection of DNA sequences. The goal of mapping in transcriptomics, epigenomics and other applications is to point out chromosomic positions either transcribed [START_REF] Sultan | A Global View of Gene Activity and Alternative Splicing by Deep Sequencing of the Human Transcriptome[END_REF], bound to a protein [START_REF] Johnson | Genome-Wide Mapping of in Vivo Protein-DNA Interactions[END_REF], or whose three dimensional conformation is altered by a protein [START_REF] Boyle | High-Resolution Mapping and Characterization of Open Chromatin across the Genome[END_REF]. Hence, further analysis only consider those reads that mapped to a unique genomic positions. In other frameworks, all mapped reads inclusive those mapped at multiple positions provide important information to detect, e.g. new copies of repeats in the sampled genome. The number of (uniquely and/or multi-) mapped reads depends on the read length, on the expected probability for a read to map on the genome, on the level of sequence errors in the reads, as well as on the genetic differences between the cell from which the reads were sequenced and that which provided the reference genome sequence. Two approaches are possible: to map exactly or approximately (up to a limited differences number between the read and the genome sequences) reads on the genome sequence. The choice between the two is not obvious since it has been shown for instance that exact mapping with a shorter read length can yield the same number of uniquely mapped reads than approximate matching up to two mismatches [START_REF] Rivals | Mpscan: fast localisation of multiple reads in genomes[END_REF], and since all approximate mapping tools are not based on the same algorithm [START_REF] Jiang | Seqmap: mapping massive amount of oligonucleotides to the genome[END_REF][START_REF] Smith | Using quality scores and longer reads improves accuracy of solexa read mapping[END_REF][START_REF] Li | SOAP: short oligonucleotide alignment program[END_REF][START_REF] Hao Lin | ZOOM! Zillions of oligos mapped[END_REF][START_REF] Langmead | Ultrafast and memory-efficient alignment of short dna sequences to the human genome[END_REF]. If approximate mapping is used, another question is then how to distinguish a difference due to a genetic variation or to sequence error in a match? More practically, whether or not sequence quality information is provided aside the reads themselves, often the complete read sequence cannot be exploited because of low quality positions. Hence, either preprocessing with various parameters is applied to eliminate some positions or multiple mapping with different parameters are tested to optimize the mapping output. In any case, the number of reads to map is so large that mapping efficiency and scalability, both in terms of time and to a less extent of memory, becomes a major issue. In Section 1.4, we will discuss on the comparison of exact versus approximate mapping approaches on these issues. Before that, Section 1.2 presents efficient solutions for the single pattern matching problem, while Section 1.3 details fast algorithms for multiple or set pattern matching.

SINGLE PATTERN MATCHING ALGORITHMS

We consider exact string matching for locating the occurrences of single nucleotide or amino acid sequence patterns in long biological sequences. We assume that the sequences are in the raw format.

Algorithms for DNA Sequences

Most of the efficient string matching algorithms in the DNA alphabet are modifications of the Boyer-Moore algorithm [START_REF] Boyer | A fast string searching algorithm[END_REF], which applies two heuristics: match and occurrence. Most often only the occurrence heuristic (also called the bad character heuristic) is applied for shifting. The Boyer-Moore-Horspool algorithm [START_REF] Horspool | Practical fast searching in strings[END_REF] (BMH) is the most famous implementation of this simplification. Because the DNA alphabet contains only four symbols, shifts based on one character are short on average. Therefore, it is advantageous to apply q-mers (or q-grams), strings of q characters, for shifting instead of single characters. This technique was already mentioned in the original paper of Boyer and Moore [11, p. 772], and Knuth [34, p. 341] analyzed theoretically its gain. Zhu and Takaoka [START_REF] Zhu | On improving the average case of the Boyer-Moore string matching algorithm[END_REF] presented the first algorithm utilizing the idea. Their algorithm uses two characters for indexing a two dimensional array.

Baeza-Yates [START_REF] Baeza-Yates | Improved string searching[END_REF] introduced an extension of the BMH algorithm, where the shift array is indexed with an integer formed from a q-mer with shift and add instructions. For this kind of approach the practical upper limit with 8 bit characters is two characters.

For the DNA alphabet Kim and Shawe-Taylor [START_REF] Kim | Fast string matching using an n-gram algorithm[END_REF] introduced a convenient alphabet compression by masking the three lowest bits of ASCII characters. In addition to the a, c, g, and t one gets distinguishable codes also for n and u. Even the important control code \n=LF has a distinct value, but \r=CR gets the same code as u. With this method they were able to use q-mers of up to six characters. Indexing of the shift array is similar to Baeza-Yates' algorithm.

With the DNA alphabet the probability of an arbitrary short q-mer appearing in a long pattern is high. This restricts the average shift length. Kim and Shawe-Taylor [START_REF] Kim | Fast string matching using an n-gram algorithm[END_REF] introduced a variation for the cases where the q-mer in the text occurs in the pattern. Then two additional characters are checked one by one to achieve a longer shift.

In most cases the q-mer that is taken from the text does not match with the last q-mer of the pattern, and the pattern can be shifted forward. For efficiency one can use a skip loop [START_REF] Hume | Fast string searching[END_REF], where the pattern is moved forward until the last q-mer of the pattern matches with a q-mer in the text. The easiest way to implement this idea is to place a copy of the pattern as a stopper after the text and artificially define the shift of the last q-mer of the pattern to be zero. After a skip loop the pattern is compared with the corresponding text positions.

A crucial thing for the efficiency of a q-mer algorithm is how q-mers are computed. Tarhio and Peltola [START_REF] Tarhio | String matching in the DNA alphabet[END_REF] presented a q-mer variation of BMH that applies a skip loop. The algorithm computes an integer called fingerprint from a q-mer. The ASCII codes are mapped to the range of 4: 0 ≤ r[x] ≤ 3, where r[x] is the new code of x, such that characters a, c, g, and t get different codes and other possible characters get e.g. code 0. In this way the computation is limited to the effective alphabet of four characters. The fingerprint is simply a reversed number of base 4. A separate transformation table h i is used for each position i of a q-mer and multiplications are incorporated during preprocessing into the tables:

h i [x] = r[x]• 4 i . For q = 4, the fingerprint of x 0 • • • x 3 is 3 i=0 r[x i ]• 4 i , which is then computed as h 0 [x 0 ] + h 1 [x 1 ] + h 2 [x 2 ] + h 3 [x 3 ].
Recently Lecroq [START_REF] Lecroq | Fast exact string matching algorithms[END_REF] presented a related algorithm. Its implementation is based on the Wu-Manber algorithm [START_REF] Wu | A fast algorithm for multi-pattern searching[END_REF] for multiple string matching, but as suggested above, the idea is older [START_REF] Boyer | A fast string searching algorithm[END_REF][START_REF] Zhu | On improving the average case of the Boyer-Moore string matching algorithm[END_REF]. For q = 4, the fingerprint of

x 0 • • • x 3 is ((((((x 0 << 1) + x 1 ) << 1) + x 2 ) << 1) + x 3 ) mod 255.
SSABS [START_REF] Sheik | A FAST pattern matching algorithm[END_REF] and TVSBS [START_REF] Thathoo | TVSBS: A fast exact pattern matching algorithm for biological sequences[END_REF] were developed with biological sequences in mind. SSABS is a Boyer-Moore type algorithm. In the search phase the algorithm verifies that the first and last characters of the pattern match with the current alignment before checking the rest of the alignment (a.k.a. guard tests). TVSBS uses a 2-mer for calculating the shift, adopted from the Berry-Ravindran algorithm [START_REF] Berry | A fast string matching algorithm and experimental results[END_REF], which is a cross of the Zhu-Takaoka algorithm and Sunday's QS algorithm [START_REF] Sunday | A very fast substring search algorithm[END_REF]. Instead of the two-dimensional shift table of Berry-Ravindran, TVSBS uses a hash function to compute an index to a onedimensional table. According to Kalsi et al. [START_REF] Kalsi | Comparison of exact string matching algorithms for biological sequences[END_REF], SSABS and TVSBS are not competitive with q-mer algorithms in the DNA alphabet.

We present one fast q-mer algorithm for DNA sequences in detail. It is SBNDM4 [START_REF] Durian | Tuning BNDM with q-grams[END_REF], which is a tuned version of BNDM (Backward Nondeterministic DAWG Matching) by Navarro and Raffinot [START_REF] Navarro | Fast and flexible string matching by combining bit-parallelism and suffix automata[END_REF]. BNDM is a kind of cross of the Backward DAWG Matching algorithm (BDM) [START_REF] Crochemore | Text algorithms[END_REF] and Shift-Or [START_REF] Baeza-Yates | A new approach to text searching[END_REF] algorithms. The idea of BNDM is similar as in BDM, while instead of building a deterministic automaton, a nondeterministic automaton is simulated even without constructing it. The resulting code applies bit-parallelism and it is efficient and compact. We present a C code for SBNDM4 as Alg. 1.1. There P[0],...,P[m-1] is the pattern and T[0],...,T[n-1] the text. The code contains a skip loop so that a copy of the pattern is placed to

Algorithm 1.1 01 for (i = 0; i < 256; i++) B[i] = 0; 02 for (i = 0; i < m; i++) B[P[m-i]] |= (1<<i); 03 for (i = 0; i < m; i++) T[n+i] = P[i]; 04 while (1) { 05 while (!(d=(B[T[j]]<<3)&(B[T[j-1]]<<2)&(B[T[j-2]]<<1)&B[T[j-3]])) 06 j +=
T[n],...,T[n+m-1]
. This version prints only the number of matches. The matches can be reported by changing line 12.

SBNDM4 is very fast in practice. One x86 processors one can still boost its performance by using 16-bit reading [START_REF] Durian | Tuning BNDM with q-grams[END_REF]. In searching DNA patterns of 20 characters, SBNDM4 with 16-bit reading is more than eight times faster than the classical Boyer-Moore algorithm [START_REF] Durian | Tuning BNDM with q-grams[END_REF] (see also comparisons [START_REF] Kalsi | Comparison of exact string matching algorithms for biological sequences[END_REF][START_REF] Lecroq | Fast exact string matching algorithms[END_REF][START_REF] Tarhio | String matching in the DNA alphabet[END_REF]).

SBNDM4 works for DNA patterns of up to 32 or 64 characters depending on the word size of the processor. In practice longer exact patterns are seldom interesting, but e.g. Lecroq's algorithm is good for them.

SBNDM4, like other variations of BNDM, works also for more general string matching where positions in the pattern or in the text represent character classes [START_REF] Navarro | Flexible Pattern Matching in Strings -Practical on-line search algorithms for texts and biological sequences[END_REF] instead of single characters. So e.g. the standard IUB/IUPAC nucleic acid codes can be used with SBNDM4.

There are also algorithms [START_REF] Rautio | String matching with stopper encoding and code splitting In[END_REF][START_REF] Kim | Fast matching method for DNA sequences[END_REF] for packed DNA. We decided to leave them outside this presentation.

Algorithms for Amino Acids

In general, there is not much difference in searching amino acid and natural language patterns. So any good search algorithm for natural language is also applicable to amino acids. In searching short patterns [START_REF] Durian | Tuning BNDM with q-grams[END_REF], the 2-mer variation of SBNMD4 is among the best. SBNDM2 is got from SBNDM4 as follows: replace line 5 by

while (!(d=(B[T[j]]<<1)&B[T[j-1]]))

and m -3 by m -1, j -4 by j -2, and m -4 by m -2 on lines 6-9.

As in the case of SBNDM4 one can still boost the performance of SBNDM2 by using 16-bit reading [START_REF] Durian | Tuning BNDM with q-grams[END_REF]. In searching patterns of 5 characters, SBNDM2

ALGORITHMS FOR MULTIPLE PATTERNS

In this section, we consider exact searching of multiple patterns. More precisely, we are given a text and r patterns and we need to find all occurrences of all the patterns in the text.

Trie-Based Algorithms

Many algorithms for exact searching of multiple patterns are based on a data structure called trie for storing the patterns. A trie is a tree where each edge is labeled with a character. Each node of the trie is associated with a string that is formed by concatenating all the labels of the edges on the path from the root to the node. Given a node in the trie, all edges to the children of this node have a different label. Figure 1.1 shows the trie storing the strings {acc, accgt, attg, ccat, cgtt }.

Aho-Corasick

The Aho-Corasick algorithm [START_REF] Aho | Efficient string matching: An aid to bibliographic search[END_REF] builds as preprocessing an automaton that recognizes the occurrences of all patterns. The preprocessing starts by building the trie of the pattern set. We add an edge from the root to the root for all those characters that do not yet have an outgoing edge from the root. The trie is then augmented with failure links as follows. The failure link of a node N in the trie points to a node that is associated with the longest possible suffix of the string associated with the node N excluding the node N itself. Additionally we associate an output function with each node whose associated string is one of the patterns. The output function outputs the identifier of the pattern. The automaton is used for searching the text as follows. We start at the root node of the trie. We read the text character by character and for each character we perform the following actions. While there is no child node with an edge labeled with the read character, we follow the failure link. Then we descend to the child with an edge labeled with the read character. Finally we output the identifiers returned by the output function for the child node.

Set Backward Oracle Matching

The Set Backward Oracle Matching (SBOM) algorithm [START_REF] Allauzen | Factor oracle of a set of words[END_REF] builds an automaton that recognizes at least all factors of the reversed patterns. The automaton is built as follows. First we build a trie of the reversed patterns. Then we traverse the trie in breath-first order and add some more edges between the nodes turning the trie into a directed acyclic graph (DAG) that recognizes all factors of the reversed patterns. To sbom preprocess (P 1 , . . . , P r ) 1. build trie(P r 1 , . . . , P r r ) 2. set supply link of root to NULL 3. for each node N in the trie in breath first order 4. down = supply link of parent 5. c = the edge label from parent to N 6.

while (down = NULL and down does not have a child with edge label c) 7.

add an edge from down to N with label c 8. down = supply link of down 9.

if (down = NULL) 10.

set supply link of N to the child of down with edge label c 11. else 12.

set supply link of N to root assist us in adding these new edges we associate a supply link with each node.

For each node we then perform the pseudo-code on lines 4 to 12 shown in Figure 1.3. Figure 1.4 shows an example of the SBOM automaton built for patterns {aag, gac} As we can see the automaton recognizes also some other strings, like caa, than the factors of the patterns. This automaton is then used for searching the occurrences of the patterns as follows. Initially we set the endpoint to the length of the shortest pattern. We then read the characters of the text backward starting at the endpoint character. For each character we make the corresponding state transition in the automaton. Whenever we encounter a node associated with one of the patterns, we verify the read region character by character against the pattern. If we encounter a character that does not have a transition on that character, we can shift the endpoint forward and start the backward scan again at the new endpoint. The length of the shift is m -j, where m is the length of the shortest pattern and j is the number of characters that we have read, or the length of the shortest pattern, whichever is shorter.

Filtering Algorithms

Filtration aims at eliminating most positions that cannot match any pattern with an easy criterion. Then, verification checks whether the remaining positions truly match a pattern. Thus filtering algorithms operate in three phases. The patterns are first preprocessed, in the second phase we search the text with a filtering method, and the candidate matches produced by the filtering are verified in the third phase.

Here we describe several algorithms that use a generalized pattern of character classes for filtration [START_REF] Salmela | Multipattern string matching with qgrams[END_REF]. Let us explain the filtration scheme with an example. Assume a set of 3 patterns of length m = 8: {P 1 , P 2 , P 3 } = {accttggc, gtcttggc, accttcca}, and set q to 5. The overlapping 5-mers (or 5grams) of each pattern are given in Figure 1.5. For a text window W of length 8 to match P 1 , the substring of length q starting at position i in W must match the i th q-mer of P 1 for all possible i, and conversely. Now, we want to filter out windows that do not match any pattern. If the substring starting at position i in W does not match the i th q-mer of neither P 1 , P 2 , nor P 3 , then we are sure W cannot match any of the patterns. Thus, our filtration criterion to surely eliminate any non-matching window W is to find if there exists a position i such that the previous condition is true.

Given a set of patterns, the filtering algorithms build a single q-mer generalized pattern (Fig. 1.5c). A generalised pattern allows several symbols to match at a position (like a position [DEN Q] in a PROSITE pattern, which matches the symbols D, E, N, and Q). However, here each q-mer is processed as a single symbol. Then, a string matching algorithm that can handle classes of characters is used for searching for occurrences of the generalized pattern in the text.

Various different algorithms can be used for implementing the filtering phase. Below we describe in more detail algorithms where filtering is based on the shift-or, BNDM, and Boyer-Moore-Horspool algorithms. A filtering algorithm always requires an exact algorithm to verify the candidate matches. In principle, any of the presented exact algorithms could be used for this purpose.

Multi-Pattern Shift-Or with q-Grams

The shift-or algorithm is easily extended to handle classes of characters in the pattern [START_REF] Abrahamson | Generalized string matching[END_REF][START_REF] Baeza-Yates | A new approach to text searching[END_REF], and thus developing a filtering algorithm for multiple pattern matching is straightforward. The {P 1 , P 2 , P 3 } = {accttggc, gtcttggc, accttcca} (a) preprocessing phase now initializes the bit vectors for each q-mer as follows.

The i:th bit is set to 0 if the given q-mer is included in the character class in the i:th position. Otherwise the bit is set to 1. The filtering phase proceeds then exactly like the matching phase of the shift-or algorithm. Given this scheme, it is clear that all actual occurrences of the patterns in the text are candidates. However, there are also false positives as the generalized pattern matches also other strings than the original patterns.

Multi-Pattern BNDM with q-Grams

The second filtering algorithm is based on the BNDM algorithm by Navarro and Raffinot [START_REF] Navarro | Fast and flexible string matching by combining bit-parallelism and suffix automata[END_REF]. This algorithm has been extended to classes of characters in the same way as the shift-or algorithm. We call the resulting multiple pattern filtering algorithm BG (short for BNDM with q-Grams). The bit vectors of the BNDM algorithm are initialized in the preprocessing phase so that the i:th bit is 1 if the corresponding q-mer is included in the character class of the reversed generalized pattern in position i. In the filtering phase, the matching is then done with these bit vectors. As with SOG, all match candidates reported by this algorithm must be verified.

Multi-Pattern

Horspool with q-Grams The last of our algorithms uses a Boyer-Moore-Horspool [START_REF] Horspool | Practical fast searching in strings[END_REF] type method for matching the generalized pattern against the text. Strictly speaking, this algorithm does not handle character classes properly. It will return all those positions where the generalized pattern matches and also some others. This algorithm is called HG (short for Horspool with q-Grams). hg matcher (T = t 1 . . . t n , n)

5-mer

1. i = 1 2. while(i ≤ n -m + 1) 3. j = m -q + 1 4.
while The preprocessing phase of HG constructs a bit table for each of the m-q+1 positions where a q-mer starts in the pattern. The first table keeps track of q-mers contained in the character class of the first position of the generalized pattern, the second table keeps track of q-mers contained in the character classes of the first and the second position in the generalized pattern, and so on. Finally, the m -q + 1:st table keeps track of characters contained in any of the character classes of the generalized pattern. Figure 1.6a shows the four tables corresponding to the pattern 'acccttaa' when using 5-mers.

These tables can then be used in the filtering phase as follows. First, the m -q + 1 st q-mer is compared with the m -q + 1 st table. If the q-mer does not appear in this table, the q-mer cannot be contained in the character classes of positions 1 . . . m -q + 1 in the generalized pattern, and a shift of m -q + 1 characters can be made. If the character is found in this table, the m -q:th character is compared to the m -q:th table. A shift of m -q characters can be made if the character does not appear in this table and therefore not in any character class in the generalized pattern in positions 1, . . . , m -q. This process is continued until the algorithm has advanced to the first table and found a match candidate there. The pseudo code is shown in Figure 1.6b. Given this procedure, it is clear that all positions matching the generalized pattern are found. However, also other strings will be reported as matches.

Other Algorithms

Other algorithms for searching multiple patterns include the Commentz-Walter algorithm [START_REF] Commentz-Walter | A string matching algorithm fast on the average[END_REF] with its variations, the Wu-Manber algorithm [START_REF] Wu | A fast algorithm for multi-pattern searching[END_REF], and algorithms derived from the Rabin-Karp algorithm [START_REF] Karp | Efficient randomized pattern-matching algorithms[END_REF] for a single pattern.

APPLICATION OF EXACT SET PATTERN MATCHING FOR READ MAPPING AND COMPARISON WITH MAPPING TOOLS

Here, we concentrate on the question of set pattern matching and on its main current application: read mapping on genomic sequences. In most frameworks, millions of reads that originate from a genome have been sequenced using HTS. A read serves as a signature for a molecule or a chromosomic position. The goal of mapping is to find back for each different read the chromosomic position of origin in the reference genome. As a read may be sequenced several time according to its number of occurrence in the biological sample, the number of different reads may be much lower than the number of read sequences. For example, in a transcriptomic assay were 2 million reads were sequenced, the read set contains ≃ 440.000 elements [START_REF] Philippe | Using reads to annotate the genome: influence of length, background distribution, and sequence errors on prediction capacity[END_REF]. As a read sequence can differ from the original chromosomic sequence because of polymorphisms or sequence errors, read mapping is often performed using approximate pattern matching, which allows for a few mismatches and or indels. For approximate mapping, either near-exact sequence similarity search programs ((BLAT [START_REF] Kent | BLAT-The BLAST-Like Alignment Tool[END_REF], MegaBlast [START_REF] Zhang | A greedy algorithm for aligning DNA sequences[END_REF], or SSAHA [START_REF] Ning | SSAHA: A Fast Search Method for large DNA Databases[END_REF]) or mapping tools (Eland, tagger [START_REF] Iseli | Indexing Strategies for Rapid Searches of Short Words in Genome Sequences[END_REF], rmap [START_REF] Smith | Using quality scores and longer reads improves accuracy of solexa read mapping[END_REF], seqmap [START_REF] Jiang | Seqmap: mapping massive amount of oligonucleotides to the genome[END_REF], SOAP [START_REF] Li | SOAP: short oligonucleotide alignment program[END_REF], MAQ [START_REF] Li | Mapping short DNA sequencing reads and calling variants using mapping quality scores[END_REF], Bowtie [START_REF] Langmead | Ultrafast and memory-efficient alignment of short dna sequences to the human genome[END_REF],

and ZOOM [START_REF] Hao Lin | ZOOM! Zillions of oligos mapped[END_REF]) are used. An alternative option when dealing with short reads is to resort to exact set pattern matching, for which mpscan offers an efficient solution [START_REF] Rivals | Mpscan: fast localisation of multiple reads in genomes[END_REF][START_REF] Philippe | Using reads to annotate the genome: influence of length, background distribution, and sequence errors on prediction capacity[END_REF].

Due the number of reads to match, repeated application of a single pattern matching algorithm for each read would require an unaffordable computing time. Hence, practically efficient solutions involve 1. either indexing the reads in main memory and scanning the genome only once (or a few times) for all reads (the solution chosen in MegaBlast [START_REF] Zhang | A greedy algorithm for aligning DNA sequences[END_REF], seqmap [START_REF] Jiang | Seqmap: mapping massive amount of oligonucleotides to the genome[END_REF], or mpscan [START_REF] Rivals | Mpscan: fast localisation of multiple reads in genomes[END_REF])

2. or first preprocessing the genome to build an index and then loading the index in memory before searching each read one after the other (the approach followed in SSAHA [START_REF] Ning | SSAHA: A Fast Search Method for large DNA Databases[END_REF], BLAT [START_REF] Kent | BLAT-The BLAST-Like Alignment Tool[END_REF], and in mapping tools like Eland, tagger [START_REF] Iseli | Indexing Strategies for Rapid Searches of Short Words in Genome Sequences[END_REF], rmap [START_REF] Smith | Using quality scores and longer reads improves accuracy of solexa read mapping[END_REF], SOAP [START_REF] Li | SOAP: short oligonucleotide alignment program[END_REF], MAQ [START_REF] Li | Mapping short DNA sequencing reads and calling variants using mapping quality scores[END_REF], Bowtie [START_REF] Langmead | Ultrafast and memory-efficient alignment of short dna sequences to the human genome[END_REF]).

1.4.1 mpscan: an efficient exact set pattern matching tool for DNA/RNA sequences

The program mpscan [START_REF] Rivals | Mpscan: fast localisation of multiple reads in genomes[END_REF][START_REF] Philippe | Using reads to annotate the genome: influence of length, background distribution, and sequence errors on prediction capacity[END_REF] is an implementation of the exact Multi-Pattern BNDM with q-Grams algorithm (cf. Section 1.3.2.2, also called BG algorithm in [START_REF] Salmela | Multipattern string matching with qgrams[END_REF]). It is specialised for searching large sets of relatively short DNA/RNA patterns in large DNA/RNA sequence files, and its interface is adequate for the purpose of mapping reads: it handles file formats commonly used in biology, can search for the reverse complementary of the pattern, etc. Its correctness, which ensures it to yield for each read all text positions at which that read matches the text, derives from that of the Multi-Pattern BNDM with q-Grams algorithm [START_REF] Salmela | Multipattern string matching with qgrams[END_REF]. The filtration efficiency depends on the parameter q. We have shown recently that the average time complexity of mpscan for searching r patterns of size l in a text of length n over an alphabet of size c is O(n log c (rl)/l) provided that q = Θ(log c (rl)).

As it was proved that the minimum time required is Ω(n log c (rl)/l) [START_REF] Navarro | Average complexity of exact and approximate multiple string matching[END_REF], this makes mpscan asymptotically optimal in average.

For example on an Intel Xeon CPU 5140 processor at 2.33 GHz with 8 GB main memory, when searching 4 million 27 bp reads on the 247 Mbp of human chromosome 1, mpscan sets the parameter q to 13, uses 229 megabytes memory, and takes 78 seconds.

Other solutions for mapping reads

With HTS becoming more popular and the increase of their sequencing capacity, the question of mapping reads on a genome sequence is a crucial issue, as well as a bottleneck.

At the HTS advent, an available solution was to use ultrafast similarity search BLAST-like programs, which were not designed for this purpose, but for locally aligning sequences that differ little (for instance, only because of sequencing errors). There were typically intended to align Expressed Sequence Tags on the human genome. These programs are not adapted to short reads (below 60 bp), and because of internal limitation cannot handle millions of queries. Hence, both their sensitivity and scalability are insufficient for mapping application with short reads [START_REF] Rivals | Mpscan: fast localisation of multiple reads in genomes[END_REF]. However, some user still resort to these tools for they allow, unlike mapping tools, an unrestricted number of differences between the read and the genome [START_REF] Peter V Kharchenko | Design and analysis of ChIP-seq experiments for DNA-binding proteins[END_REF]. All these tools implement a filtration strategy that requires a substring of the query sequence to match the genome either exactly [START_REF] Ning | SSAHA: A Fast Search Method for large DNA Databases[END_REF][START_REF] Zhang | A greedy algorithm for aligning DNA sequences[END_REF] or with at most one mismatch [START_REF] Kent | BLAT-The BLAST-Like Alignment Tool[END_REF].

Since the commercialisation of HTS, numbers of commercial or free mapping tools have been developped or published (cf. list above); for instance the Eland software is provided with the Illumina R Solexa sequencer. As mentioned the goal of mapping differs with the application but it is often to find the best match for a read: the match with the least differences, and if possible unique. All mapping programs perform successive approximate pattern matching up to a limited number of differences. Some tools can find matches with up to 4 mismatches and/or indels, but generally a guarantee to find all matches (as required in the definition of approximate matching) is given only up to one or two mismatches. This limitation makes sense to speed up the search and derives from the applied filtration scheme. All tools (except ZOOM) use variant of the so-called PEX filter [START_REF] Navarro | Flexible Pattern Matching in Strings -Practical on-line search algorithms for texts and biological sequences[END_REF], which consists in splitting the read in k + 1 adjacent pieces, knowing that at least one piece will match exactly when a maximum of k errors are allowed. Many mapping programs makes it efficient by using 2-bit encoded sequences and/or an index of the genome (e.g., MAQ, Eland, Bowtie, rmap, SOAP).

The program ZOOM exploits spaced seeds: it requires that a subsequence of a defined form, instead of a substring, matches between the read and the genome [START_REF] Burkhardt | Better filtering with gapped q-grams[END_REF][START_REF] Ma | PatternHunter: faster and more sensitive homology search[END_REF]. The subsequence's pattern of required matching positions and wild-cards is designed on purpose depending on the expected match length and maximal number of differences [START_REF] Kucherov | Multiseed Lossless Filtration[END_REF]. The advantage of spaced seeds is their capacity to handle mismatches and insertion/deletions (indels), and their increase sensitivity compared to substring based filtration [START_REF] Burkhardt | Better filtering with gapped q-grams[END_REF][START_REF] Ma | PatternHunter: faster and more sensitive homology search[END_REF]. Their main drawback is the difficulty of seed design: ZOOM uses a conjunction of several seeds. Hence, sets of spaced seeds are specifically designed for a certain read/match length and a maximum number of allowed differences, and different sets corresponding to different parameter combinations are hard coded in ZOOM. All known formulations of the seed design problem are at least NP-hard, even for a single seed [START_REF] Kucherov | Multiseed Lossless Filtration[END_REF][START_REF] Ma | On the complexity of the spaced seeds[END_REF][START_REF] Nicolas | Hardness of optimal spaced seed design[END_REF].

Comparison of mapping solutions

As already mentioned, many groups have developed and/or published their own mapping tools, and all tools, except mpscan, implement a solution based on approximate pattern matching. However to date, one lacks a comparative evaluation of the sensibility of all these tools in various application frameworks. The intended application makes a difference since for e.g., identifying genomic variations, multiple matching locations of a read provides useful information, while in transcriptomics one usually discards multi-mapped reads. Probing the sensitivity and evaluating the sensitivity vs speed or memory balance is a difficult task knowing that the programs differ in their notion of approximation (e.g., with or without indels).

Here, we discuss the conclusions of a comparison on the less difficult task of exact set pattern matching. We exclude the program Eland for it is not free for academics, as well as MAQ, which does not accept parameters for searching only exact read matches. 1.4.3.1 Speed, memory footprint, and scalability We compared rmap, seqmap, SOAP (v1 & v2), ZOOM, Bowtie and mpscan for searching increasing read sets on the longest human chromosome (chromosome 1, 247 Mbp). The public input data sets contains 6.5 million of 27 bp reads and we took subsets every million reads (available on the GEO database under accession number GSM325934). At the date of this comparison, this set belongs to the largest ones in terms of number of different reads, and there is no available data set of similar size with much larger reads (say > 36).

Figure 1.7 reports the running times in seconds on a logarithmic scale for searching the subsets of 1, 2, . . . up to 6 and 6.5 million reads. Of course the times do not include the index construction time for those programs that use an index, which in the case of, e.g. Bowtie, lasts hours for the complete human genome. First, all tools can handle very large read sets and their running times remain impressive even if they degrade somehow with increasing read sets. Second, the comparison of ZOOM or mpscan compared to genome indexing tools like Bowtie or SOAP shows that high performances are not bound to a genome index, at least for exact pattern matching. Knowing that ZOOM algorithm also handles approximate matches with up to two mismatches or indels, it seems that it offers a very satisfying solution compared to Bowtie, which is limited to mismatches and offers less guarantee. For exact pattern matching, the performance differences can be quite large (times 10 between mpscan and SOAP-v2 for two million reads), and mpscan offers the fastest solution overall, even if it exploits only a 32-bit architecture. However, mpscan time increases more when going from 4 to 5 million reads, suggesting that for equal read length, a coarse grain parallelization would improve its performances.

To illustrate the low memory footprint of mapping tools that do not load a genome index in RAM, we give the amount of RAM required by ZOOM, seqmap, and mpscan for searching the complete human genome with one million 27 bp tags. ZOOM requires 17 minutes and 0.9 Gigabytes, rmap takes 30 min and 0.6 Gb, seqmap performs the task in 14 min with 9 Gb, while mpscan needs < 5 min using 0.3 Gb. In contrast, Bowtie human genome index, which is implemented as a Burrows-Wheeler Transform, takes at least 1.4 Gb [START_REF] Langmead | Ultrafast and memory-efficient alignment of short dna sequences to the human genome[END_REF].

1.4.3.2 Exact pattern matching for read mapping The read length influences the probability of a read to map on the genome, and also its probability to map once. The shorter the read the higher the probability of mapping, but the lower that of mapping once. In many applications, reads mapping at unique genomic positions are preferred. A rationale for the currently developped extension of read length is the increase probability to map a unique genomic location. On the human genome, a length of 19 bp already brings the risk of mapping at random below 1% and we have shown recently that it already maximise the number of uniquely mapped reads on four real data sets [START_REF] Philippe | Using reads to annotate the genome: influence of length, background distribution, and sequence errors on prediction capacity[END_REF]. Studying the sequence error position in the reads, we could show that the error probability at one sequence position increases with the position in the read for Illumina R /Solexa data. Hence, an alternative to approximate mapping is to perform exact matching using only a prefix of each read (of an adequate length).

To evaluate this we compared the result of approximate matching with full length reads with that of mpscan on read prefixes. Eland searches the best read match up to two mismatches, while we ran mpscan to search for exact matches of read prefixes. The full length read are 34 bp. If one maps with mpscan the full length reads, 86% remain unmapped and 11% are uniquely mapped. With at most two mismatches, Eland finds 14% of additional uniquely mapped reads with one or two mismatches, while mapping the 20 bp prefix of each read with mpscan allows to map 25% of all reads at unique positions (14% more sites than with full length reads). Hence, both approaches yield similar output, but exact matches represent easier and more secure information than approximate matches. For the current rates of sequencing errors and read lengths, exact matching is a suitable solution for read mapping. Moreover, it allows us to estimate computationally the sequence error rate without performing control experiments (cf. [START_REF] Philippe | Using reads to annotate the genome: influence of length, background distribution, and sequence errors on prediction capacity[END_REF] for a more in-depth presentation), which would be more difficult using approximate matching.

Conclusions

About pattern matching, q-gram based algorithms and especially mpscan represent the most efficient theoretical and practical solutions to exact set pattern matching for huge pattern sets (above a million patterns). Compared to known solutions surveyed seven years ago in [START_REF] Navarro | Flexible Pattern Matching in Strings -Practical on-line search algorithms for texts and biological sequences[END_REF], which were reported to handle several hundred thousands patterns, mpscan provides more than an order of magnitude improvement: it allows to process at astonishing speed pattern sets of several millions reads. The second take home message is that its filtration scheme can compete with approaches that use a text index.

Since 2005, the capacity of HTS is continuously evolving: biotechnological research and development aim at reducing the quantity of biological extract, augmenting the sequencing capacity and quality, raising the read length, and even enlarging the application fields. Despite the efforts for designing scalable and efficient mapping programs, it will remain a computational issue to let mapping solutions fit the requirements of new HTS versions. This is a complex question since read length above 20 are not necessary to point out a unique position in a genome as large as that of human [START_REF] Philippe | Using reads to annotate the genome: influence of length, background distribution, and sequence errors on prediction capacity[END_REF].

An interesting conclusion is that different filtration schemes achieve impressive efficiency and scalability, but may be insufficient for tomorrow's needs. The abundant pattern matching literature may still contain other possible algorithms whose applications in this setup has not been yet evaluated. With the spread of multi-core computers, parallelization represents another future line of research.

Finally, we left aside the problem of mapping pairs of reads. In this framework, two reads are sequenced for each targeted molecule: each at a different extremity. The reads come in pairs and the goal of mapping is to find one matching position for each read such that the two positions are on the same chromosome and in an upper bounded vicinity. In other applications, the pair relations are unknown and it is then required to find across the two sets of beginning and ending reads, which ones constitute a pair for they map on the same chromosome not too far from another [START_REF] Rivals | Transcriptome Annotation using Tandem SAGE Tags[END_REF]. Some mapping tools like MAQ or ZOOM can solve read pair mapping efficiently, while a precursor of mpscan has been developped and applied in the second framework [START_REF] Rivals | Transcriptome Annotation using Tandem SAGE Tags[END_REF].

Figure 1 . 1 :

 11 Figure 1.1: An example trie storing the strings {acc, accgt, attg, ccat, cgtt }

Figure 1 . 2 :

 12 Figure 1.2: An example of the Aho-Corasick automaton for the patterns {acc, accg, att, cca, cgt }. The dashed lines show the failure links. Failure links to the starting state of the automaton have been omitted.

  Figure 1.2 shows an example of an Aho-Corasick automaton for the patterns {acc, accg, att, cca, cgt }.

Figure 1 . 3 :Figure 1 . 4 :

 1314 Figure 1.3: Preprocessing of SBOM

Figure 1 . 5 :

 15 Figure 1.5: (a) A set of 3 patterns of length m = 8. (b) The overlapping 5-mers starting at position 1 to 4 (resp. in magenta, orange, red, blue) of each pattern. (c) The generalised 5-mer pattern for the set of tags.

Figure 1 . 6 :

 16 Figure 1.6: The HG algorithm: (a) the data structures for the pattern "acccttaa" and (b) the pseudo code for the search phase.

Figure 1 . 7 :

 17 Figure 1.7: Comparison of mapping tools: Search times of rmap, seqmap, SOAP (v1 & v2), ZOOM, Bowtie and mpscan in seconds (log scale) for increasing subsets of 27 bp reads. All tools behave similarly and offer acceptable scalability. mpscan remains the most efficient of all, and can be 10 times faster than tools like seqmap or rmap. Times do not include the index construction time.

A base pair is the length unit of a DNA/RNA sequence.