Florence Tubach

Dominique Salmon

Philippe Ravaud

Yannick Allanore

Philippe Goupille

Maxime Bréban

Béatrice Pallot-Prades

Sophie Pouplin

Antoinette Sacchi

Rose Marie Chichemanian

Risk of tuberculosis is

Keywords: Adult, Aged, Antibodies, Monoclonal, adverse effects, therapeutic use, Arthritis, Rheumatoid, drug therapy, Behcet Syndrome, drug therapy, Case-Control Studies, Colitis, drug therapy, Female, France, Humans, Immunoglobulin G, adverse effects, therapeutic use, Male, Middle Aged, Prospective Studies, Receptors, Tumor Necrosis Factor, therapeutic use, Registries, Risk Factors, Spondylarthritis, drug therapy, Treatment Outcome, Tuberculosis, chemically induced, drug therapy, epidemiology, Tumor Necrosis Factor-alpha, immunology anti-TNF-alpha, tuberculosis, safety, rheumatoid arthritis, inflammatory chronic diseases

à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Risk of tuberculosis is higher with anti-tumor necrosis factor monoclonal antibody therapy than with soluble tumor necrosis factor receptor therapy: The three-year prospective French Research Axed on Tolerance of Biotherapies registry

Florence Tubach 2 1 # , Dominique Salmon 3 # , Philippe Ravaud [START_REF] Askling | Risk and case characteristics of tuberculosis in rheumatoid arthritis associated with tumor necrosis factor antagonists in Sweden[END_REF]

Service d h pato-gastro-ent rologie 12

Methods

An incidence study with the French population as reference and a case-control analysis. We collected, for 3 years, cases of TB among French patients receiving anti-TNF therapy, whatever the indication, with two controls treated with anti-TNF agents per case.

Results

We collected 69 cases of TB in patients treated for rheumatoid arthritis (n 40), spondylarthropathies (n 18), inflammatory colitis (n All three available TNF antagonists have been associated with increased incidence of TB. However, a possible difference between infliximab and etanercept in incidence of TB is suspected from the spontaneous reporting of TB cases to the US Food and Drug Administration ( ) before September 2002, when no recommendations of screening for latent TB existed. In one study, the incidence of 5 TB was 28 per 100,000 patient-years for patients receiving etanercept versus 54 per 100,000 for patients receiving infliximab.( ) 9

= = = 9
Likewise, four registries found small differences in incidence of TB between RA patients receiving infliximab or etanercept, although not statistically significant, perhaps because in each of these registries, the number of TB cases was low (2 to 13), these registries being not powered enough to demonstrate very rare events or investigate a difference in risk depending on anti-TNF agent.( , , , ) No specific 1 2 6 7 epidemiologic data are available for adalimumab, but results of phase-3 and -4 clinical studies suggested that use of this drug may increase the risk of TB to a level close to that of infliximab ( ). 10

Thus, we used another strategy to investigate incidence of TB by setting up a national registry to collect all cases of TB occurring during 3 years in French patients receiving anti-TNF therapy, whatever the indication for use. The aim of this large prospective study was to describe cases of TB associated with anti-TNF therapy and their outcome; identify the risk factors of TB occurrence in patients receiving anti-TNF therapy, particularly the anti-TNF agent; and estimate the global and drug-specific incidence of TB after the establishment and widespread knowledge by clinicians of recommendations for TB prophylaxis before anti-TNF therapy.

PATIENTS AND METHODS

The French RATIO (Research Axed on Tolerance of bIOtherapies) registry was designed by a multidisciplinary group to prospectively collect all cases of TB occurring from February 1, 2004, to January 1, 2007, in patients who were receiving anti-TNF therapy.( ) This 11 registry involved clinicians from all concerned medical specialties, the French agency of drugs Agence Fran aise de S curit Sanitaire des

ç é é
Produits de Sant (AFSSAPS), and its network of 31 regional pharmaco-vigilance centres.

é

We conducted an incidence study and a case-control analysis to investigate the risk of newly diagnosed TB associated with the use of anti-TNF agents.

Setting

French guidelines ( ) elicited in 2002, recommended anti-TB chemoprophylaxis before all anti-TNF therapy for patients presenting 12 at least one of the following: tuberculin skin test (TST) diameter 10 mm (revised to 5 mm in July 2005), abnormal chest X-ray results

≥ ≥

with calcifications > 1 cm, previous untreated exposure to TB or episode of TB. Correct chemoprophylaxis was defined as a 9-month treatment with isoniazide or a 3-month treatment with 2 anti tuberculosis drugs including rifampicine ( ). 12

Identification and validation of tuberculosis cases

To enhance exhaustiveness of the collection of cases of TB, different sources were used.

All cases reported to the 31 French pharmaco-vigilance regional centres were collected (reporting of severe drug-related adverse events is mandatory for clinicians in France) but also cases pointed out directly to the companies commercializing anti-TNF.

In addition, physicians from all the different French hospital centres implied either in the prescription of TNF blockers (i.e. rheumatology, internal medicine, gastroenterology and dermatology departments) and/or in the management of TB (i.e. infectious diseases centres, intensive care units, chest medicine units) were directly required to report each newly diagnosed case. A direct mail reminder 4 times a year and several communications at congresses or in specialized press encouraged them to report cases.

We obtained also cases from the Institut National de Veille Sanitaire (InVS), the French national public health surveillance system that monitors TB in France, with a mandatory clinician-and biologist-based notification system.

A detailed standardized case report form was completed for each reported case of TB (whatever the source of notification), and additional documents were collected (hospitalization summary, histological and microbiological results, and other data if necessary) to enter the validation process.

Validation of cases

Included in the RATIO registry were all cases (from all sources) with a validated diagnosis of TB according to the 10 revision of the th International Classification of Diseases (A15 A19). An expert committee involving 3 experts of TB (SB, DS, OL) validated cases by consensus on the basis of the detailed standardized case report form and the additional documents.

Assessment of exposure to anti-TNF therapy

Infliximab has been available since 1999, etanercept since 2002 and adalimumab since 2004. For all, it is mandatory in France that the first prescription be from a hospital physician. For patients who received more than one anti-TNF agent, we considered three ways of addressing the exposure: first anti-TNF agent, last anti-TNF agent, and any use of the anti-TNF agent.

Incidence study

Incidence rate of TB

We estimated the annual incidence rate of TB in patients treated with anti-TNF therapy, adjusted for age and sex, with the French population used as reference. The numerator consisted in the validated cases of TB from the RATIO registry. For the denominator, we estimated the number of patient-years of receipt of anti-TNF agents in France during the 3-year period of the study (2004, 2005 and 2006) from different sources. For this purpose, we obtained the total number of doses of each anti-TNF agent sold in France for the 3 years study, both according to the national French agency of drugs (AFSSAPS), and to the 3 pharmaceutical firms (Abbott, Schering-Plough and Wyeth). Each firm provided its estimation of the number of patient-years for the 3 anti-TNF agents on the basis of the number of doses sold and the mean dosage used in each indication. These data were very consistent. We used the mean of these different estimates to derive the estimation of the number of patients receiving each anti-TNF agent that we used in the main analysis.

In addition, for sensitivity analysis, we obtained independent data from the R gime Social des Ind pendants (RSI), the French health é é insurance fund for self-employed workers, which involves 3 million enrolees (i.e., 5 of the French population representative of the whole % population) and provides claims and demographic data for use of etanercept and adalimumab.

Standardized Incidence Ratios

To assess the distribution of underlying diseases, sex and age of treated patients, we used surveys of samples of prescribing clinicians, hospital pharmacies, and nonhospital pharmacies performed by the firms, and the RSI files for the patients receiving etanercept and adalimumab.

On the basis of these data and the French annual incidence rate of TB by 10-year age and sex class from 2002 to 2006 (from the Institut National de Veille Sanitaire) ( ), we estimated the standardized incidence ratio ( ) to compare the TB risk associated with 13 14

anti-TNF treatment with that in the French population.

Risk of TB in patients receiving anti-TNF therapy

A case-control study was performed with cases from the RATIO registry.

Cases

Cases were all validated cases of TB in the registry with a labelling indication for use of anti-TNF treatment (i.e. RA, spondylarthropathy SpA; AS or psoriatic arthritis , UC or CD, or psoriasis).

[ ]

Controls

Tuberculosis-free patients receiving anti-TNF treatment in a labelling indication were included from centres involved in the RATIO registry in a global pool of controls. From that pool, we randomly selected patients for a database of controls reflecting the proportion of patients receiving each of the three anti-TNF drug in France. Two controls per case were randomly matched by sex and underlying inflammatory disease from this database of controls. We also used a second sample of controls randomly selected from the same database of controls, with the same matching criteria (second matching). In addition, we used a third sample of controls (1 control per case), adding the time from onset of anti-TNF treatment as matching criteria (third matching).

Statistical analysis

The number of cases of TB in France during the study period determined the sample size. A descriptive analysis was performed for the whole sample.

Annual incidence rate of tuberculosis

The SIR was calculated for anti-TNF agent use as a whole and for agents used individually. The 95 confidence intervals (CIs) for % each SIR were calculated, and the deviation of the SIR from 1.00 was tested assuming a Poisson distribution for the observed cases because of the small number of cases.

Case-control study

We identified the risk factors for TB by both univariate and multivariate analysis (conditional logistic regression model indicating TB sequellae at screening, last anti-TNF agent received, time from onset of anti-TNF treatment, use of a potentially immunosuppressive disease-modifying anti-rheumatic drug during the last year, use of steroids during the last year, and use of > 10 mg/day or boluses of steroids during the last year.

Factors included in the stepwise multivariate model were those associated with the status (case or control) in univariate analysis with a significance level of p<0.20. We used the bootstrap method to examine the stability of our final model (1000 replicated bootstrap samples to find candidate variables which were more often retained). A p<0.05 was considered statistically significant.

Subgroup analyses

We performed subgroup analyses for the more frequent underlying inflammatory diseases: RA and SpA

Sensitivity analyses

To better assess the drug-specific risk, we performed sensitivity analyses for the SIR and the OR measuring the risk of adalimumab or infliximab versus etanercept, varying the definition of exposure (all anti-TNF agents, by anti-TNF agent, and for patients receiving only one anti-TNF agent, receiving anti-TNF therapy during at least 6 weeks, and currently receiving anti-TNF therapy), the different denominator estimates (denominator estimates from each of the 3 pharmaceutical firms and from the RSI), and the different criteria used for matching the controls with cases (second and third matching).

The main analyses were independently performed by two statisticians. Statistical analysis involved use of SAS release 9.1 (SAS Inst., Cary, NC).

Compliance with research ethics standards

This study authorized by the ethic committee of AP-HP, GHU Nord (Institutional Review Board of Paris North Hospitals, Paris 7

University, AP-HP; authorization number 162-08). The registry was reported at clinicaltrials.gov (ClinicalTrials.gov Identifier: NCT00224562).

RESULTS

Description of the TB cases

We collected data on 75 cases of TB, and 69 cases were validated. Characteristics of the 69 patients with a validated TB diagnosis are in . In addition, one patient had diabetes, 3 a past history of cancer, none was HIV positive, and 1 had a chronic obstructive

Outcome

Three patients required hospitalization in an intensive care unit. After a median follow-up of 22.9 months, 2 patients died (one with disseminated TB died before diagnosis; another died from a cause unrelated to TB), 45 were considered cured, 5 were still receiving treatment, and in 5 TB relapsed (after a median delay of 9.4 months) (outcome was missing for 12 cases). None of the 5 patients had restarted anti-TNF therapy before relapse.

At TB diagnosis, anti-TNF therapy was not stopped in 2 patients (treated with infliximab). Moreover, anti-TNF therapy was restarted in 8 patients (etanercept in 3, infliximab in 5) after a median anti-TB treatment duration of 5.8 months (range 5 days to 10 months). No recurrence of TB was observed in these 10 patients (median follow-up 17.7 months, range 8.2 to 32.2 months).

Time occurrence of TB with anti-TNF therapy

The median time to occurrence of TB since start of anti-TNF treatment was 12.0 months. The median time since start of the last anti-TNF was 9.9 months. In 4 patients, TB had occurred despite discontinuation of anti-TNF therapy 7.1 to 13.9 months before. In these 4 patients, the anti-TNF was infliximab in 4 and adalimumab in 1.

As indicated in

, the incidence of TB was higher during the first year of anti-TNF treatment and differed depending on type of Figure 1 treatment: monoclonal antibodies (infliximab or adalimumab) or soluble receptor (etanercept).

Incidence and risk of TB for patients receiving anti-TNF therapy versus the general population

The main analysis relied on a total number of 57,711 patient-years of use of anti-TNF therapy during the 2004 2006 period, with 18 -% receiving adalimumab, 51 etanercept and 31 infliximab as the denominator of the incidence rate. The annual incidence rate of TB % % adjusted for age and sex for patients receiving anti-TNF therapy, with the French population as a reference, was 116.7 (95 CI 10.6 222.9

% per 100,000 patient-years). The SIR was 12.2 (9.7 15.5; p<0.0001).

-

The risk of TB for patients receiving anti-TNF therapy as compared to the French population differed depending on anti-TNF agent used. The annual adjusted-incidence rate of TB was 9.3 (0.0 29.4) per 100,000 for patients receiving etanercept, 187. for adalimumab. Because the incidence rates and SIRs were close for adalimumab and infliximab and their mechanisms of action are similar, we pooled data for adalimumab and infliximab for subgroup and sensitivity analyses. Still, a difference between monoclonal-antibody and soluble-receptor anti-TNF therapy was observed in the main analysis and in the different sensitivity analyses (

), even when we separately used the different estimates from independent sources, which gave very consistent adjusted incidence Figure 2 rates and SIRs ( ). suppl Figure 1

Risk factors of TB for patients receiving anti-TNF therapy

The case-control study involved 68 cases (as described in the methods, the patient treated for Behcet s disease was not included in this ' case-control study) and 136 controls (from a pool of controls whose anti-TNF treatments reflected that of all users nationally, i.e. 18% receiving adalimumab, 51 etanercept and 31 infliximab). The results of the univariate case-control analysis are given in . On

% %

Table 2 multivariate analysis, the following factors were predictive of TB: age, being born in an endemic area, time from onset of anti-TNF treatment, and therapy with adalimumab or infliximab versus etanercept.

Results of the main multivariate analysis are presented in . Risk of TB was higher with infliximab or adalimumab than with table 3 etanercept. (OR 13.3 (95 CI 2.6 69.0; p 0.002) and 17.1 (95 CI 3.6 80.6; p 0.0003), respectively. As for the SIRs, we pooled data for

= % - = % - =
adalimumab and infliximab for the subgroup and sensitivity analyses, for more power. The drug-specific OR were also similar with the second and third matching and in subgroup and sensitivity analysis ( ). Figure 3 DISCUSSION Our analysis of 69 cases of TB collected over 3 years in the RATIO registry clearly demonstrated a difference in risk of TB between patients receiving monoclonal-antibody than soluble-receptor anti-TNF therapy. We observed a 7-to 17-fold difference in risk depending on the 2 types of anti-TNF agents used, from two different methodological approaches to assess the risk of TB: a comparison of incidence rate with the French general population and a case-control study with controls receiving anti-TNF therapy for the same underlying inflammatory diseases. The risk of TB was higher during the first year of anti-TNF treatment, which favoured the reactivation of latent TB. Moreover, none of the 69 cases had been treated with correct chemoprophylaxis against TB before anti-TNF treatment was initiated.

After the FDA alert of the risk of TB associated with the use of infliximab, numerous countries set up registries to investigate the safety of anti-TNF agents. All registries of biologics but ours are cohort studies involving only a part of the population of focus and thus are not powered enough to demonstrate very rare events or investigate a difference in risk depending on anti-TNF agent used, nor are the randomized controlled trials.( ) In the cohort of the US national data bank for rheumatic diseases (NDB), the rate of TB was not 15 increased in anti-TNF-na ve RA patients (6.2 (1.6 34.4) per 100,000 patient-years) but was increased in RA patients receiving infliximab ï -(52.5 (14.3 134.4) per 100,000 patient-years)( ). However in this study, the rate was calculated from only 4 cases of TB among 6,460

- 7
patients receiving infliximab. In the same study, no case of TB was observed for 2,327 patients receiving etanercept. In the British registry, the incidence of TB was 1.5 per 1000 patient-years with infliximab and 0.5 per 1000 patient-years with etanercept, but the number of cases was 7 and 2, respectively.( ) In the Swedish registry, the incidence rates were 1.5 per 1000 patient-years with infliximab and 0.8 per 1000 2 patient-years with etanercept, 9 cases treated with infliximab alone and 4 with etanercept alone.( ) Finally, among Korean patients, 2 1 cases of TB were observed among 90 patients receiving infliximab and no case was observed among 103 patients receiving etanercept.( ) 6

The difference of risk between the types of anti-TNF agents was only suggested in these registries and is clearly demonstrated in the RATIO registry. In our registry, patients were included when the adverse event occurred, not when the therapy began. Thus, our population of focus was the whole French population, and this allowed us to collect in 3 years many more cases of TB than that collected by all the other registries. Furthermore, as far as we know, ours is the only biologic registry collecting safety data for patients receiving anti-TNF therapy, for any indication.

Contrary to the standard pharmacovigilance notifications, all the files of TB patients included in the RATIO registry were reviewed and validated by an expert committee, which is a main quality issue, especially for the diagnosis of TB, which may be difficult.

As always, incidence rates might be somewhat under-or overestimated because of an imprecise numerator or denominator.

Concerning the numerator, passive surveillance data may underestimate the incidence of adverse events, and thus the annual rate of TB in the RATIO registry might be toward the low end of the true incidence rate.( , ) The denominator of the incidence rate was estimated 16 17 only. The pharmaceutical firms may have used different methodologies for evaluating the denominator. However, because each firm evaluated the number of patient-years during the period for each anti-TNF agent, the difference in risk between agents we observed cannot be explained by the different methodologies used for the different agents. Furthermore, in the sensitivity analysis, we used the different estimates from independent sources separately, which gave very consistent adjusted incidence rates and SIRs (

). In supplementary Figure 1 addition, our estimates for adalimumab and etanercept did not differ from data obtained from the RSI, which has claims data for outpatients receiving etanercept or adalimumab among 3 million RSI enrolees --5 of the French population. The consistency of results % of sensitivity analyses from these independent sources is a strong argument for their validity. Furthermore, the difference in risk between agents we found (7-to 17-fold difference) is unlikely to be explained by a lack of precision in the denominator estimate.

Although this national survey design has some limitations (especially the lack of precision in the denominator), it is probably the only way (or at least the most powerful way) to investigate difference in risk with use of anti-TNF agents.

The mechanism by which TNF antagonists reactivate latent TB is not fully understood. In animal models, TNF plays a central role in Two-thirds of our cases of TB occurred in patients with negative TST results at screening. However, if we strictly follow the national recommendation for prophylactic treatment in France, 56.5 of the patients with a TST result of less than 5 mm and 39.5 with a result % %

of less than 10 mm should have received an anti-TB prophylactic treatment. Concerning the others, either the screening was not adequate or new TB infection developed, not reactivating latent TB. Our study suggests that taking into account whether patients were born in an endemic area could improve the screening of latent TB. However this information may have been less precisely reported in the control group. Other ways to improve the screening of latent TB, such as use of specific blood tests assessing the presence of specific in vitro anti-TB T-cells (Quantiferon TB or T-spot-TB), should be evaluated.

Our study also confirmed that, as emphasized by others,( ) anti-TB chemoprophylactic treatment is efficient for preventing TB. None 8 of our patients with TB had received a correct prophylactic antibiotic against TB. The efficacy of screening and anti-TB chemoprophylaxis may explain why the median time from the beginning of anti-TNF therapy to TB onset is much higher (12 and 10 months for the first and the last anti-TNF, respectively) in our study than in the first study describing TB in patients treated with infliximab (14 weeks),( ) when 9 no screening for TB was recommended. However, probably we observed both features of TB infection in our study: reactivation of latent TB infection (corresponding to the high risk during the first year of anti-TNF treatment as shown in and in the case-control Figure 1 results) and new infection occurring during follow-up, which will never be avoided with pre-treatment screening.

In conclusion, our study provides clear evidence of a higher risk of TB with monoclonal-antibody anti-TNF therapy than with soluble-receptor therapy. These differences can be supported by differences of action of the 2 types of agents on membrane TNF leading to a differential effect on effector T cells and on Treg. These different mechanisms of action could also explain a better efficacy of monoclonal-antibody therapy in CD, in other granulomatous diseases such as sarcoidosis( ) and in uveitis.( ) 28 29

Figure 1

Cumulative incidence of TB as a function of anti-TNF treatment duration as a whole and by anti-TNF agent.

Figure 2

Analysis of the estimation of the SIR for the risk of TB, whatever the anti-TNF agent, and according to the anti-TNF agent received and the underlying inflammatory disease n is the number of cases involved in the calculation (numerator of the incidence rate).The plot size relates to the number of patients treated involved in the calculation (denominator of the incidence rate).

Figure 3

Sensitivity analysis of the results of the case-control analysis odds ratios (OR) from multivariate analysis for the risk of being treated with adalimumab or infliximab rather than with etanercept. n is the number of cases involved in the calculation. The plot size relates to the number of patients treated involved in the calculation.

¥ Time from onset of last anti-TNF treatment and first symptoms of TB for the cases, time from onset of last anti-TNF treatment and last news for controls * Immunosuppressive disease-modifying anti-rheumatic drugs (DMARDs): at least one among cellcept, chlorambucil, ciclosporin, cychlophosphamid, mercaptopurine, mycophenolate mofetil, methotrexate, leflunomide, azathioprine, abatacept, anakinra, cortolizumab, CTLA4 Ig, rituximab. 0.001 OR: odds ratio ¥ Time from onset of last anti-TNF treatment and first symptoms of tuberculosis for cases; time from onset of last anti-TNF treatment and last news for controls

INTRODUCTION

  ), psoriasis (n 1) and Beh et s disease (n 1) treated with infliximab (n 36), adalimumab (n 28) and etanercept (n 5). None of the = ç ' = = = = cases had received correct chemoprophylaxis treatment. The sex and age-adjusted incidence rate of TB was 116.7 per 100,000 patient-years. The SIR was 12.2 (95 confidence interval 9.7 15.5) and was higher for therapy with infliximab and adalimumab than % -Treatment with tumor necrosis factor (TNF) antagonists has been recognized as a risk factor of active tuberculosis (TB) in patients with immune-mediated inflammatory diseases, including rheumatoid arthritis (RA), ankylosing spondylitis (AS), Crohn s disease (CD), ' psoriatic arthritis, and psoriasis.( ) Most cases of TB develop soon after treatment initiation and correspond to reactivation of a latent 1 -7 TB infection.( , ) Scientific organizations and health authorities worldwide have proposed recommendations for screening patients with 4 5 latent TB infection and treating such patients before initiating anti-TNF treatment, and the effectiveness of these recommendations has been demonstrated.( ) 8

  mediating mycobacterial infections, and soluble TNF but especially membrane TNF, expressed by activated macrophages and T lymphocytes, is essential in protecting against TB infection.( ) In a murine model of chronic TB, neutralization of TNF by antibody but 18 not the TNF-receptor fusion molecule exacerbates chronic TB because of better penetration of antibodies into the granuloma.( ) Humans show higher avidity and better stability in membrane TNF with anti-TNF monoclonal antibodies,( ) which leads in some studies to more 20 efficient apoptosis.( ) Even if in other studies, infliximab and etanercept showed no clear difference in inducing apoptosis ( 21 -23 in vivo ) the difference in reverse signalling due to a difference in membrane TNF targeting may have other consequences. For example, 24 anti-TB-specific T cells show a clear difference in effector function. Two studies provided the same results: adjunct infliximab or in vitro adalimumab treatment at therapeutic concentrations to specific human anti-TB T cells was much more efficient in decreasing proliferation and interferon gamma secretion by these T cells than adjunct etanercept treatment at therapeutic concentrations.( , ) As well, 20 25 inhibition of interleukin 1 released by lipopolysaccharide-stimulated monocytes is greater with monoclonal-antibody anti-TNF treatment than with soluble-receptor treatment.( ) Finally, monoclonal anti-TNF therapy alleviates Treg impairment seen in RA by inducing 26 apoptosis of membrane TNF-positive nonefficient Treg, thus allowing re-expansion of efficient regulatory T cells (Treg).( ) Thus, our 27 epidemiological results can be supported by differences in action of the two types of anti-TNF agents on specific effector T-cells and on Treg.

2 , Yannick Allanore 4 , Philippe Goupille 5 , Maxime Br ban é 6 , Bé atrice Pallot-Prades 7 , Sophie Pouplin 8 , Antoinette Sacchi 9 , Rose Marie Chichemanian 10 , St phane Bretagne é 11 , Dominique Emilie 12 , Marc Lemann 13 , Olivier Lorthololary 14 , Xavier Mariette 15 * for the RATIO group

	1	Mod les et m thodes de l valuation th rapeutique des maladies chroniques è é 'é é	INSERM : U738 , Universit Paris-Diderot -Paris VII é	, Facult
	é	de m decine Paris 7 16, Rue Henri Huchard 75018 Paris,FR é
	D partement d pid miologie, biostatistique et recherche clinique 2 é 'é é	AP-HP , H pital Bichat -Claude Bernard ô	, 46 rue Henri Huchard 75018
	Paris,FR
	3	Service de m decine interne é	AP-HP , H pital Cochin ô	, Universit Paris Descartes é	, 27 rue du Faubourg Saint Antoine, 75679 Paris Cedex
	14, France,FR
	4	Service de rhumatologie	AP-HP , H pital Cochin ô	, Universit Paris Descartes é	, 27 rue du Faubourg Saint-Jacques 75679 Paris Cedex
	14,FR
	5	Service de rhumatologie	CHU Tours , Universit Fran ois Rabelais -Tours é ç	, Tours,FR
	6	Service de Rhumatologie		AP-HP , H pital Ambroise Par ô	é , 9 avenue Charles-de-Gaulle 92100 Boulogne,FR
	7	Service de rhumatologie	H pital Bellevue

ô , Saint-Etienne,FR Service de rhumatologie 8 CHU Rouen , 76000 Rouen,FR CH F. Quesnay, Service de m decine interne 9 é Centre hospitalier F. Quesnay , Mantes la Jolie,FR Centre r gional de pharmacovigilance 10 é Centre r gional de pharmacovigilance é , Nice,FR Service de parasitologie 11 AP-HP , H pital Henri Mondor ô , Universit Paris-Est é , Cr teil,FR é

Table 1 pulmonary

 1 disease. For the 58 patients receiving only one anti-TNF agent, 34 (58.6 ) had received infliximab, 23 (39.7 ) adalimumab Risk factors of TB, before anti-TNF therapy, were available in 51 patients. 34 (66.7 ) showed at least one risk factor (TST diameter % ≥ 5 mm, chest x-ray results showing signs of a history of TB, history of an untreated TB or of exposure to TB, or being born in an endemic area): the TST diameter, available in 45 patients, was > 10 mm for 4 patients (9 ), 5 to 10 mm for 11(24 ), < 5 mm for 30 (67 ). For 6 of 57 patients (10.5 ) (not done in 12), the chest x-ray findings suggested a history of TB. In addition, 10 patients had a history of % exposure to TB, 1 a personal history of TB, and 3 a history of primary TB infection, and 14 were born in an endemic area. According to French guidelines for TB chemoprophylaxis before anti-TNF initiation, with a TST diameter cut-off of 5 or 10 mm, 26 of 46 (56.5 ) or 17 % of 43 (39.5 ) of the patients should have received chemoprophylaxis, respectively. None of the patients had received correct

	%	%	%

% % and 1 (1.7 ) etanercept. Two cases (and 4 controls) had also received anakinra.

%

Table 3

 3 Risk factors of TB for patients receiving anti-TNF therapy (multivariate analysis)

	TB cases and controls (68 cases and 136 controls)
	OR 95 CI [ %	]	p value

[START_REF] Kelly | Failed validation of a clinical decision rule for the use of radiography in acute ankle injury[END_REF]

[START_REF] Kelly | Failed validation of a clinical decision rule for the use of radiography in acute ankle injury[END_REF]

[START_REF] Kelly | Failed validation of a clinical decision rule for the use of radiography in acute ankle injury[END_REF]

[START_REF] Kelly | Failed validation of a clinical decision rule for the use of radiography in acute ankle injury[END_REF]

[START_REF] Kelly | Failed validation of a clinical decision rule for the use of radiography in acute ankle injury[END_REF]

[START_REF] Kelly | Failed validation of a clinical decision rule for the use of radiography in acute ankle injury[END_REF]

[START_REF] Kelly | Failed validation of a clinical decision rule for the use of radiography in acute ankle injury[END_REF]

[START_REF] Kelly | Failed validation of a clinical decision rule for the use of radiography in acute ankle injury[END_REF]

[START_REF] Kelly | Failed validation of a clinical decision rule for the use of radiography in acute ankle injury[END_REF]

Ackowledgements:

The authors thank all the clinicians who had reported at least one case in the RATIO Registry: Abitbol (Paris), Allanore (Paris), André (Clermont-Ferrand), Ardizzone (Mulhouse), Bergman (Paris), Azais (Poitiers), Bachelez (Paris), Bardet (Orl ans), Beau (Poitiers), Bergman é (Paris), Belmatoug (Clichy), Berthelot (Nantes), Blasco (Barbois), Bonnet (Limoges), Bouhnik (Clichy), Bourgarit (Paris), Bouvard (Angers), Bressot (Chalon sur Saone), Brian on (Aix les bains), Brocq (Nice), Cadiot (Reims), Castela (Nice), Visanica (Metz), Combert (La Rochelle), ç Couret (Valence), Cuillerier (Dreux), Dalle (Lyon), Dasilva (Elbeuf), Debandt (Aulnay sous-bois), Debourdeau (Lyon), Depernet (Langers), Dereure (Montpellier), Descamps (Paris), Duriez (Saint-Brieuc), Fach (Bergerac), Fain (Bondy), Fautrel (Paris), Filippi (Nice), Fior (Bondy), Flourie (Lyon), Fulpin (Marseille), Gaborit (Orange), Gaudin (Grenoble), Gendre (Paris), Ghringhelli (Bordeaux), Gillet (Nancy), Goupille (Tours), Grados (Amiens), Grosclaude (Uriage), Gueit (Rouen), Guillaume (Colmar), Guyot (Roubaix), Hersebach (Rennes), Hoen (Besanç on), Houvennagel (Lomme), Beguinot (Reims), Jang-Guyro (Brian on), Jardin (Privas), Justum (Caen), Laharie (Bordeaux), Lambotte (Le ç kremlin-Bic tre), Lecompte (Nancy), Leparc (Boulogne-Billancourt), Lequen (Pau), L hirondel (Caen), Lin (Les Soissons), Liot (Paris), The authors are also grateful to G.R. Auleley, J Deligne and C Blum-Boisgard from the RSI for providing data to validate the denominator estimate of the incidence rate we used, D Che and D Antoine from the Institut National de Veille Sanitaire for providing the French annual incidence rate of tuberculosis by age and sex, and the AFSSAPS and the regional pharmacovigilance centres for their contribution to the exhaustiveness of the RATIO registry.

The authors thank C Roy and G Baron for the statistical analysis, and N Nicolas, S Makhlouf and A Djemoui for their help in collecting and preparing the validation of the cases.

Role of the funding source

The pharmaceutical companies (Abbott, Schering Plough and Wyeth) had no role in the study design, data collection, data analysis, data interpretation or writing of the report. The corresponding author had full access to all the data and had final responsibility for the decision to submit for publication.

Allanore received consulting and/or talk honorarium from Wyeth (less than 10,000 ), Goupille received consulting and/or talk honorarium

±

Categorized variables are in no, ( ).

% § Time from onset of first/last anti-TNF treatment to first symptoms of TB, time from onset of last anti-TNF treatment.

Tuberculosis complicating anti-TNF therapy Arthritis Rheum . Author manuscript Page / 13 14 § Time from onset of first/last anti-TNF treatment to first symptoms of TB, time from onset of last anti-TNF treatment. * immunosuppressive disease-modifying anti-rheumatic drugs (DMARDs): at least one among cellcept, chlorambucil, ciclosporin, cychlophosphamid, mercaptopurine, mycophenolate mofetil, methotrexate, leflunomide, azathioprine, abatacept, anakinra, cortolizumab, CTLA4 Ig, rituximab. Continuous variables are in mean (SD). Categorized variables are in no, ( ).

% OR: odds ratio ¥ Time from onset of last anti-TNF treatment and first symptoms of TB for the cases, time from onset of last anti-TNF treatment and last news for controls


domain: spi math phys

text:
Rémy Charayron

Thierry Lefebvre
email: thierry.lefebvre@onera.fr

Nathalie Bartoli
email: nathalie.bartoli@onera.fr

Joseph Morlier
email: joseph.morlier@isae-supaero.fr

Multi-fidelity Bayesian optimization strategy applied to Overall Drone Design

Nowadays, drones can be developed for a wide range of use cases, from infrastructure monitoring to sea rescue, urban mobility or military purposes. Which drone design is best suited for a specific mission? To answer this question, we need to solve a constrained optimization problem based on a multi-disciplinary design model that takes the mission into account. The model generally being a computationally expensive numerical model whose gradients are not available all the time encourages us to consider a Bayesian optimization approach. Such strategy is well known to achieve a trade-off between exploitation and exploration in order to find interesting minimal area with a reduced number of function evaluations. A multi-fidelity approach can improve even more the computational efficiency of the Bayesian optimization strategy. In this work, we aim at designing a fixed-wing drone (fully electric) for long range surveillance mission. Two fidelity level electric drone models are developed. For a given mission requirement, the final battery state of charge is optimized with respect to drone design variables. Optimizations are performed on several missions using both a mono and a multi-fidelity Bayesian optimization strategy. The interest of using a multi-fidelity method for overall drone design has been assessed. The multi-fidelity super-efficient global optimization algorithm (MFSEGO) appeared to need less budget to reach convergence than the mono-fidelity algorithm and to be more robust to the initial design of experiments.

I. Nomenclature

II. Introduction

O ptimization of systems with strong multidisciplinary interactions is one of the current challenges in the aircraft design community. This context of multidisciplinary design analysis and optimization (MDAO) [START_REF] Martins | Engineering Design Optimization[END_REF] is reflected in this work. Of the various disciplines that could fit into a multidisciplinary Unmanned Aerial Vehicle (UAV) model, the mission discipline is important but not the easiest to consider at a detailed level. In this paper, we have endeavored to achieve a high fidelity mission discipline. Regarding the optimization strategy, Bayesian technique [START_REF] Frazier | A tutorial on Bayesian optimization[END_REF] was chosen to solve the problem, it is a well-known approach in the field of expensive black box optimization as it does not need the gradient information. In this work, we focus on the value of multi-fidelity. When multiple sources of information are available, we can take advantage of the characteristics of each source. In general, cheap but inaccurate information sources are used to explore the design space, while expensive and accurate sources are used for exploitation. First, a brief state of the art on multi-fidelity Kriging and Bayesian optimization is provided in Section III. Next, Section IV describes the two developed high-fidelity (HF) and low-fidelity (LF) UAV models that take into account the mission. Section V presents the four test cases experimented in this paper. Each test case corresponds to a geometric drone design for a specific mission. In Section VI, the results for each test case are outlined. Finally, the conclusion and perspectives are given in Section VII.

III. State of the art

This section presents the multi-fidelity framework used in this paper. We want to solve an optimization problem with inequality or/and equality constraints:

min 𝑥 ∈Ω 𝑓 (𝑥) such that 𝑔 𝑖 (𝑥) ≤ 0 𝑔 𝑒 (𝑥) = 0 (1)
where the objective function 𝑓 , the inequality constraint functions 𝑔 𝑖 and the equality constraint functions 𝑔 𝑒 are evaluated with a Multidisciplinary Design Analysis (MDA) costly to evaluate. In a multi-fidelity context, the 𝑓 , 𝑔 𝑖 and 𝑔 𝑒 can be obtained by different MDAs associated to different fidelity levels. To solve this multidisciplinary optimization (MDO) problem with expensive black-box functions, Bayesian optimizers [START_REF] Frazier | A tutorial on Bayesian optimization[END_REF][START_REF] Archetti | Bayesian optimization and data science[END_REF] are used to minimize the number of function evaluations. First Section III.A shows how to construct a multi-fidelity Gaussian process surrogate model. Then, in Section III.B, the current Bayesian optimization approach is adapted to multi-fidelity by defining a criterion to choose the fidelity level to query. We assume here that the reader is familiar both with Gaussian processes (GP) interpolation also denoted Kriging models [START_REF] Rasmussen | Gaussian processes for machine learning[END_REF][START_REF] Krige | A statistical approach to some basic mine valuation problems on the Witwatersrand[END_REF][START_REF] Matheron | Mémoires du Bureau de Recherche Géologiques et Minières[END_REF] and with the classical mono-fidelity Kriging based Bayesian optimization methods like the Efficient Global Optimization (EGO) [START_REF] Jones | Efficient global optimization of expensive black-box functions[END_REF] algorithm for unconstrained problems and the Super Efficient Global Optimization (SEGO) algorithm [START_REF] Sasena | Flexibility and efficiency enhancements for constrained global design optimization with kriging approximations[END_REF][START_REF] Bartoli | An adaptive optimization strategy based on mixture of experts for wing aerodynamic design optimization[END_REF][START_REF] Bartoli | Improvement of efficient global optimization with application to aircraft wing design[END_REF] for constrained problems.

A. Multi-fidelity Kriging

The first extension of the Kriging model for multivariate functions, called co-Kriging, was first developed in geostatistics (see [START_REF] Chiles | Geostatistics: modeling spatial uncertainty[END_REF][START_REF] Wackernagel | Multivariate geostatistics: an introduction with applications[END_REF]). Next, [START_REF] Kennedy | Bayesian calibration of computer models[END_REF] figured that making assumptions to relate different levels of fidelity is a way to simplify multi-fidelity problems. In this section, focus will be made on Le Gratiet's [START_REF] Gratiet | Multi-fidelity Gaussian process regression for computer experiments[END_REF] recursive formulation that we used to construct multi-fidelity Kriging models. Let's suppose that 𝐿 + 1 fidelity levels 𝑓 0 , ..., 𝑓 𝐿 sorted from the lowest to the highest are available to approach the function 𝑓 (the same philosophy could be applied to the constraint functions 𝑔 𝑖 and 𝑔 𝑒 ) and that we dispose of 𝐿 + 1 design of experiments (DoE), one for each fidelity level, denoted (𝐷 𝑖 ) 𝑖=0...𝐿 with 𝐷 𝑖 = ((𝑥 0 , 𝑓 𝑖 (𝑥 0 )), ..., (𝑥 𝑗 , 𝑓 𝑖 (𝑥 𝑗 )), ..., (𝑥 𝑁 𝑖 , 𝑓 𝑖 (𝑥 𝑁 𝑖 ))) 𝑇 , 𝑁 𝑖 being the number of points in 𝐷 𝑖 . We assume that we build some nested DoE: each point evaluated at a fidelity level is also evaluated at all the lowest fidelity levels: 𝐷 𝑥 𝐿 ⊆ 𝐷 𝑥 𝐿-1 ⊆ . . . ⊆ 𝐷 𝑥 0 where 𝐷 𝑥 𝑖 = (𝑥 0 , ..., 𝑥 𝑗 , ..., 𝑥 𝑁 𝑖 ) 𝑇 corresponds to the inputs part of 𝐷 𝑖 . The nested DoE provides some properties useful to express the surrogate model variance in closed form. The following assumption first introduced by [START_REF] Kennedy | Bayesian calibration of computer models[END_REF] in the bi-fidelity case is used to link the fidelity levels:

𝑓 𝑙+1 (𝑥) = 𝜌 𝑙 𝑓 𝑙 (𝑥) + 𝛿 𝑙+1 (𝑥) such that 𝑓 𝑙 ⊥ 𝛿 𝑙+1 ∀𝑙 = 0, ..., 𝐿 -1 (2)
where 𝛿 𝑙+1 (𝑥) is a discrepancy function that captures the difference between the 𝑙 + 1-th and the 𝑙-th fidelity levels while 𝜌 𝑙 is a scaling factor applied to 𝑓 𝑙 . Le Gratiet [START_REF] Gratiet | Multi-fidelity Gaussian process regression for computer experiments[END_REF] proposed to add the lowest fidelity function to the basis function set (ℎ 𝑖 ) 𝑖=1... 𝑝 used in the universal Kriging regression term to get:

𝜇(𝑥) = ∑ 𝑖=1,..., 𝑝 𝛽 𝑖 ℎ 𝑖 (𝑥) + 𝛽 𝜌 0 𝑓 0 (𝑥) (3)
where 𝛽 𝜌 0 is an estimation of 𝜌 0 and (𝛽 𝑖 ) 𝑖=1... 𝑝 is a set of unknown coefficients to be multiplied by the basis functions.

Estimations of these coefficients are done by maximizing the likelihood [START_REF] Bachoc | Cross validation and maximum likelihood estimations of hyper-parameters of Gaussian processes with model misspecification[END_REF][START_REF] Pavlyuk | Computing the maximum likelihood estimates: concentrated likelihood, EM-algorithm[END_REF]. Since we used a nested DoE structure, the independence between fidelity levels of the surrogate model is assumed. Then the following recursive formulation [START_REF] Gratiet | Multi-fidelity Gaussian process regression for computer experiments[END_REF] for the mean ( μ𝑙+1 ) and variance ( σ2 𝑙+1 ) of each fidelity level GP surrogate model can be written:

∀𝑙 = 0, ..., 𝐿 -1 μ𝑙+1 = 𝜌 𝑙 μ𝑙 + μ 𝛿 𝑙+1 σ2 𝑙+1 = 𝜌 2 𝑙 σ2 𝑙 + σ2 𝛿 𝑙+1 (4)
In this case, 𝜌 𝑙 is considered as a constant but it can depend on 𝑥 and we have 𝜌 𝑙 : 𝑥 ↦ → 𝜌 𝑙 (𝑥). The learning process is the following, the lowest fidelity level is first learnt, then the relationships (scaling factor 𝜌 𝑙 and discrepancy function 𝛿 𝑙+1 with 𝑙 = 0, ..., 𝐿 -1) between each fidelity level are successively learnt. Since the variances can be expressed in a closed form, the contribution of each fidelity level to the total variance of the multi-fidelity model can also be expressed. This is the most remarkable advantage of Le Gratiet's recursive formulation. Denoting 𝜎 2 𝑐𝑜𝑛𝑡 (𝑙, 𝑥) the variance contribution of the 𝑙 𝑡 ℎ fidelity level at the point 𝑥, with the notation 𝜎 2 𝛿 0 = 𝜎 2 0 , we have:

∀𝑙 = 0, ..., 𝐿 -1 𝜎 2 𝑐𝑜𝑛𝑡 (𝑙, 𝑥) = 𝜎 2 𝛿 𝑙 (𝑥) 𝐿-1 𝑗=𝑙 𝜌 2 𝑗 and 𝜎 2 𝑐𝑜𝑛𝑡 (𝐿, 𝑥) = 𝜎 2 𝛿 𝐿 (𝑥) (5)

B. Multi-fidelity Bayesian optimization

With a multi-fidelity Bayesian optimization process, one must decide the most promising point and the fidelity level at which to evaluate it. It has been proposed in [START_REF] Gratiet | Multi-fidelity Gaussian process regression for computer experiments[END_REF] to solve the problem of finding the point and the fidelity level in two successive steps. First the point is found using a classical acquisition function as in mono-fidelity Bayesian optimization. Then the knowledge of the variance contribution to each fidelity level gives some information to smartly decide which fidelity level to choose. Let 𝑐 0 , ..., 𝑐 𝐿 be respectively the querying costs of all the fidelity levels 𝑓 0 , ..., 𝑓 𝐿 . Let us denote 𝜎 2 𝑟 𝑒𝑑 (𝑙, 𝑥 * ) the variance reduction of the high fidelity model when the point 𝑥 * is evaluated with all the fidelity levels ≤ 𝑙 to ensure the nested DoE structure

𝜎 2 𝑟 𝑒𝑑 (𝑙, 𝑥 * ) = 𝑙 ∑ 𝑖=0 𝜎 2 𝛿 𝑖 (𝑥 * ) 𝐿-1 𝑗=𝑖 𝜌 2 𝑗 (6)
A criterion to choose the level of enrichment can then be written as:

𝑙 * = arg max 𝑙 ∈0,...,𝐿 𝜎 2 𝑟 𝑒𝑑 (𝑙, 𝑥 * ) ( 𝑙 𝑖=0 𝑐 𝑖 ) 2 (7)
One can remark that Eq [START_REF] Jones | Efficient global optimization of expensive black-box functions[END_REF] is a trade off between variance reduction and computational cost. Then the cost ratios between the different fidelity levels are important parameters to the method. This two step approach from [START_REF] Charayron | Multi-fidelity constrained Bayesian optimization, application to drone design[END_REF][START_REF] Meliani | Multi-fidelity efficient global optimization: Methodology and application to airfoil shape design[END_REF] is described in Fig. 11 of Appendix ??. It is refered in the following as the Multi-Fidelity Super Efficient Global Optimization (MFSEGO). Now that we have detailed the multi-fidelity model and the optimization technique that we will use in this work, we describe the electric drone models we have developed for the purposes of this work.

IV. Electric drone models developed

The increasing use of drones and in particular electric drones for various missions is accompanied by a need to model these drones and the missions they are assigned. In this paper, we focus on a power line monitoring mission. We will consider that the mission is a succession of checkpoints along the power line through which the drone must pass. Two models of electric drones have been developed in this work: a high fidelity model and a low fidelity one. Both of them couple different disciplines (aerodynamics, structure, masses, geometry, mission, propulsion, ...) which are computed with more or less accuracy. The purpose of the models is to estimate different quantities of interest: the final state of charge of the drone's battery (𝑆𝑂𝐶) at the end of the mission as well as the mission time, and two other outputs that indicate if the wing or the tail can be broken during the mission. Both models must be able to consider the UAV geometrical parameters such as wing span, root chord, taper ratio and wing and tail sweep angle, dihedral angle and spanwise twist control points as variables. In addition, structural parameters like wing and tail spanwise spar thickness control points are also considered as variables. Finally we want to model a function 𝑓 depending on 19 design variables and described in the following equation.

𝑓 : Ω ⊂ R 19 -→ R 4 𝑥 ↦ -→ 𝑦 (8)
All the models have been developed in Python using OpenMDAO, an open source MDO framework [START_REF] Gray | OpenMDAO: An open-source framework for multidisciplinary design, analysis, and optimization[END_REF] from NASA Glenn. The Extended Design Structure matrix (XDSM) [START_REF] Martins | Engineering Design Optimization[END_REF] diagrams of our OpenMDAO analyses have been made thanks to WhatsOpt [START_REF] Lafage | WhatsOpt: a web application for multidisciplinary design analysis and optimization[END_REF] which is a web application developed at ONERA to define and share multidisciplinary analyses in terms of disciplines and data exchange. Figure 1 is an example of XDSM that represents the LF model workflow. The aerostructural discipline relies on a tool called OpenAeroStruct (OAS) [START_REF] Jasa | Open-source coupled aerostructural optimization using Python[END_REF][START_REF] Chauhan | Low-fidelity aerostructural optimization of aircraft wings with a simplified wingbox model using OpenAeroStruct[END_REF][START_REF] Chaudhuri | Multifidelity optimization under uncertainty for a tailless aircraft[END_REF]. It is an open source lightweight tool that performs aerostructural analysis and optimization using OpenMDAO. It couples a vortex-lattice method (VLM) and a 6 degree of freedom 3-dimensional spatial beam model to simulate aerodynamic and structural analyses using lifting surfaces. OpenAeroStruct is developed by the NASA and the University of Michigan. In the following, two fidelity levels are presented to model the drone. For the aerostructural discipline OpenAeroStruct is used for both fidelity levels with different meshes. To deal with the mission in the multidisciplinary process, the low fidelity model relies on the electrical Breguet range equation while the high fidelity model requires to solve an optimal control problem which is tackled by the DYMOS software [START_REF] Falck | dymos: A Python package for optimal control of multidisciplinary systems[END_REF]. DYMOS is a framework for the simulation and optimization of dynamical systems within the OpenMDAO environment. It has two primary objectives: provide a generic ordinary differential equation (ODE) integration interface that allows for the analysis of dynamical systems and allow the user to solve optimal control problems involving dynamical multidisciplinary systems. In Section IV.A the low fidelity drone model is described and in Section IV.B the high fidelity one is detailed.

A. Low fidelity electrical drone model

The low fidelity model tries to approach the final 𝑆𝑂𝐶, the mission time, and the wing and tail failure constraints with a cheap computational cost. Figure 1 depicts the associated workflow with the different disciplinary modules involved in the low fidelity model and their associated number of inputs and outputs variables. Regarding performance during mission computation, this model relies on an adaptation of the Breguet range equation to electrical aircraft. The classic fuel Breguet range equation is recalled in Appendix VII.B. The electric adaptation is given by:

𝑅 = 𝐶 𝐿 𝐶 𝐷 𝜂 𝑝 𝜂 𝑚 𝜂 𝐸𝑆𝐶 𝑒 𝑏 𝑔 𝑚 batt 𝑚 (9)
with 𝜂 𝑝 , 𝜂 𝑚 and 𝜂 𝐸𝑆𝐶 being respectively the efficiency coefficients for the propeller, the motor and the electronic speed controller. 𝑒 𝑏 is the battery specific energy, 𝑔 the gravitational acceleration, 𝑚 batt the battery mass, and 𝑚 is the total mass of the aircraft. Typically, for our study, we choose 𝑒 𝑏 = 250 Wh.kg -1 . This equation is valid in the steady level flight case and its demonstration is given in Appendix VII.C.

Fig. 1 XDSM for the low fidelity drone model with the different disciplinary modules involved.

The low fidelity workflow is therefore the following. First, the RootChordTail group (see Fig. 1) balances the tail root chord in order to ensure longitudinal stability for the drone. Figure 2 shows the XDSM of the RootChordTail group. Note that in order to ensure the longitudinal stability of the drone, we use a rule of thumb. The moment arm between the tail and the wing is equal to 3 meters. The tail root chord is balanced so that the tail surface is equal to 0.15 times the wing surface. In order to obtain a realistic tail root chord, we impose that the taper ratios of the wing and tail are equal and that the wing span is equal to 3 times the tail span. The tail span, taper ratio and root chord are not input variables as they depend on the geometric properties of the wing. Next, the aerostructure group (see Fig. 1) uses OpenAeroStruct to compute the lift and drag coefficient required for the electrical Breguet component, the wing and tail failures that will be used as constraints by the optimizer and the wing and tail masses. The material used for the wing and tail spars is carbon fiber, and the following properties have been defined in the aerostructure model: the Young modulus of the spars is equal to 85𝑒9Pa, the shear modulus of the spars is equal to 25𝑒9Pa, the maximum authorized yield stress is set to 350𝑒6 𝑚 Pa with 𝑚 = 2.5 being a safety margin and the material volumic mass is equal to 1.6𝑒3 kgm -3 . In addition, we choose a tubular spar for both the wing and the tail. The radius of these spars is considered constant in the spanwise direction with a radius of 0.044m for the wing tube and a radius of 0.02m for the tail tube. The thickness to chord ratio of the wing and of the tail are chosen to be constant in the spanwise direction and equal to 0.12. Note that this value of thickness to chord ratio ensures that the tubular spars can be effectively contained within the wing and the tail. It is also important to remark that we have chosen a fixed discretization of the wing and tail meshes in the OpenAeroStruct model. Indeed the number of points in the spanwise direction is equal to 9 for the wing and to 7 for the tail while the number of points in the chordwise direction is equal to 7 for the wing and to 5 for the tail. The wing and tail masses from the aerostructure component are passed to the masses component which deduces the overall drone mass. Finally the electrical Breguet range equation is applied to find the maximum range and the functions component (see Fig. 1) is used to make a rough approximation of the 𝑆𝑂𝐶 at the end of the mission estimated by the following formula:

𝑆𝑂𝐶 = (1 - 𝑅 target 𝑅 max
)100 [START_REF] Bartoli | Improvement of efficient global optimization with application to aircraft wing design[END_REF] where 𝑅 target is the targeted range of the mission and 𝑅 max is the maximum range that the drone can perform (computed with the electrical Breguet range equation Eq ( 9)).

B. High fidelity electrical drone model

The high fidelity model also attempts to approximate the final 𝑆𝑂𝐶 and the failure constraints for the wing and the tail, but with greater accuracy. Figure 3 shows the XDSM diagram of the high fidelity drone model with the main disciplinary modules involved.

Fig. 3 XDSM for the high fidelity drone model with the different disciplinary modules involved.

The strategy involved computes the optimal trajectory in terms of energy consumption and then returns the final 𝑆𝑂𝐶 passing through this optimal path. The workflow of the high fidelity model is as follows. First, as with the low fidelity model, the tail root chord is balanced to ensure the longitudinal stability of the UAV. Next, we perform the structure sizing using a call to OpenAeroStruct at the most critical flight point: maximum speed and maximum angle of attack. This component computes the wing and tail failure constraints and the wing and tail masses. Note that, as for the low fidelity model, the wing and tail spars are considered tubular. The spar radius, the spar material properties, the thickness to chord ratio, and mesh discretization are also considered to be the same as in the low fidelity model for the wing and the tail. Finally, the mission is completed and the corresponding component provides the mission time and the battery state of charge. The main part of the model is encapsulated in the mission component which calls DYMOS [START_REF] Falck | dymos: A Python package for optimal control of multidisciplinary systems[END_REF], the OpenMDAO package for optimal control. The mission component has 𝑛 control points of the trajectory through which the drone must pass and which define the mission requirements. Each control point is composed of 4 values: altitude, range, Mach number, and climb rate. We then build 𝑛 -1 phases between each consecutive control point. To be consistent with the low fidelity model, the trajectory control points are chosen so that the sum of each of the distances of the 𝑛 -1 phases is equal to the target range 𝑅 target of the low fidelity model. DYMOS is used to determine the detailed optimal trajectory that minimizes the final 𝑆𝑂𝐶 on each phase knowing the set of design variables. In short, DYMOS looks for the optimal control variables (climb rate and Mach number) that maximize the final 𝑆𝑂𝐶 for a fixed UAV design with the constraint of passing through the trajectory control points on each phase. A constraint is added to the optimal control problem in order to ensure that the maximum possible thrust is not exceeded during the mission. Denoting 𝑛 phases the number of phases, the optimal control problem solved with DYMOS on each phase is described by the formula: ∀𝑖 ∈ 0, . . . , 𝑛 phases , 𝑢 * 𝑖 = max

𝑢 𝑖 ∈𝑈 𝑖 ∫ 𝑡 𝑖+1 𝑡 𝑖 𝜕𝑆𝑂𝐶 𝜕𝑡 (𝑠 𝑖 (𝑡), 𝑢 𝑖 (𝑡))𝑑𝑡 such that                                Thrust(𝑡) ≤ Thrust max ∀𝑡 ∈ [𝑡 𝑖 , 𝑡 𝑖+1 ] alt(𝑡 𝑖 ) = alt 𝑖 range(𝑡 𝑖 ) = range 𝑖 Mach(𝑡 𝑖 ) = Mach 𝑖 climb rate(𝑡 𝑖 ) = climb rate 𝑖 alt(𝑡 𝑖+1 ) = alt 𝑖+1 range(𝑡 𝑖+1 ) = range 𝑖+1 Mach(𝑡 𝑖+1 ) = Mach 𝑖+1 climb rate(𝑡 𝑖+1 ) = climb rate 𝑖+1 (11)
where 𝑠 𝑖 and 𝑢 𝑖 are respectively states and controls over the 𝑖-th phase and (alt 𝑖 , range 𝑖 , Mach 𝑖 , climb rate 𝑖 ) is the 𝑖-th trajectory control point. Note that we perform a maximization since the state of charge time derivative is always negative. We consider the control system associated with our optimal control problem. It can be written as a system of differential equations of the following form:

𝑠 = 𝐻 (𝑠, 𝑢) (12)
The control variables 𝑢 are the climb rate and the Mach number and the state variables 𝑠 are the 𝑆𝑂𝐶, the altitude and the range. 𝐻 is a function of the state and control variables which computes the time derivatives 𝑠 of the state variables.

To solve the optimal control problem, the time derivatives of all the state variables must be computed. To do so DYMOS allows the time derivatives to be provided either as control variables or as outputs of an OpenMDAO model able to compute them by approximating the 𝐻 function defined by Eq. ( 12). To approximate the state variables time derivatives 𝑠, we proposed a model denoted by Ĥ (.), coded in OpenMDAO, and depicted in Fig. 4.

Fig. 4 XDSM of the model used to compute the optimal control problem state variable time derivatives. It is used in the mission component which relies on DYMOS.

The altitude time derivative is the climb rate which is given as a control variable. To compute the range time derivative, we first compute the altitude dependent characteristics of the atmosphere. Then, we compute the true air speed (𝑇 𝐴𝑆) using the Mach number and the speed of sound in the air in these atmospheric conditions. Knowing the 𝑇 𝐴𝑆 and the climb rate, we can apply elementary trigonometric considerations to compute the flight path angle 𝛾 and the range rate which is the speed component along the 𝑥-axis. The 𝑆𝑂𝐶 rate ( 𝑆𝑂𝐶) is the most complex part of the model Ĥ (.) to compute. It is computed in the propulsion group whose XDSM diagram is given in Fig. 5. The 𝑆𝑂𝐶 rate is the opposite of 𝑃 batt the power output of the battery at each instant scaled by the product between the battery capacity 𝑄 batt (Coulomb) and the battery voltage 𝑈 batt (Volt). We express 𝑆𝑂𝐶 in percentages: where 𝑈 batt = 𝑈 cell × 𝑘, with 𝑈 cell being the electric tension of a single cell in the on-board Lithium-Polymer (LiPo) battery and 𝑘 being the number of cells in the battery. Typically, 𝑈 cell = 3.7𝑉 and 𝑘 = 3. 𝑃 batt is necessarily greater than the propeller output power 𝑃 prop required for flight equilibrium due to the efficiencies of the propeller, motor and electronic speed controller (𝐸 𝑆𝐶) being less than 1. These efficiency losses passing through the propeller, motor and 𝐸 𝑆𝐶 are modeled using efficiency coefficients ∈ [0; 1[, respectively 𝜂 𝑝 = 0.7, 𝜂 𝑚 = 0.95 and 𝜂 𝐸𝑆𝐶 = 0.98 and we have:

𝑆𝑂𝐶 (𝑡) = 𝜕𝑆𝑂𝐶 𝜕𝑡 (𝑡) = -𝑃 batt 100 𝑄 batt 𝑈 batt (13)
𝑃 prop = 𝑃 batt 𝜂 𝑝 𝜂 𝑚 𝜂 𝐸𝑆𝐶 (14)
The required power output of the propeller depends on the required thrust 𝑇 for flight balance and the true air speed 𝑇 𝐴𝑆

𝑃 prop = 𝑇 × 𝑇 𝐴𝑆 (15)
with the thrust 𝑇 being equal to the product of dynamic pressure 𝑞 = 1 2 𝜌 air 𝑇 𝐴𝑆, the reference area 𝑆 = 𝑆 wing + 𝑆 tail and the thrust coefficient

𝐶 𝑇 𝑇 = 𝑞𝑆𝐶 𝑇 (16)
The dynamic pressure 𝑞 is computed in a dedicated component inside the Ĥ model while the reference surface is given to the Ĥ model as an input coming from the structure component. The 𝐶 𝑇 is an output of the steady flight equilibrium group of Ĥ. Note that the trajectories considered are simple, so not dealing with acceleration phases and only assuming steady flight phases is accurate enough for us. Within the steady flight equilibrium group, we first call the aerodynamic component. Using OpenAeroStruct simulations directly at this stage to compute the aerodynamics is very expensive in term of computational cost because a very large number of simulations must be performed. Then, for the UAV defined by the MDA inputs variables, two surrogate models of the OpenAeroStruct aerodynamic model are built to solve the optimal control problem at a reasonable cost. The same radii, material properties, the same thickness to chord ratios and same mesh discretizations as in the low fidelity model are considered. These surrogate models are Kriging models built at the very beginning of the mission component from a 90 point Latin Hypercube Sampling (LHS) design of experiments. Once built the surrogate models are passed to Ĥ through an option dictionary to be used in it. They take three inputs (altitude, Mach number and angle of attack). The first surrogate provides the lift coefficient 𝐶 𝐿 while the second provides the drag coefficient 𝐶 𝐷 .

𝐶 𝐷 = 𝐶 𝐷,𝐺 𝑃 (AoA, Mach, alt) ( 17)

𝐶 𝐿 = 𝐶 𝐿,𝐺 𝑃 (AoA, Mach, alt) (18)
Next we use the 𝐶 𝐷 to compute the 𝐶 𝑇 required for steady flight equilibrium before computing the 𝐶 𝐿 required for steady flight equilibrium. We denote by 𝑊 total the total weight of the UAV. It is computed in a previous component of Ĥ the time derivative model and is equal to the sum of the battery weight (fixed equals to 15 kg), the payload weight (fixed equals to 10 kg), the fuselage weight (fixed equals to 20 kg) and the lifting surface structural weight (computed by the structural component of the model).

𝐶 𝑇 = 𝑊 total sin(𝛾) cos( 𝐴𝑜 𝐴)𝑞𝑆 + 𝐶 𝐷 cos( 𝐴𝑜 𝐴) ( 19
)
𝐶 𝐿 𝑒𝑞 = 𝑊 total cos(𝛾) 𝑞𝑆 -𝐶 𝑇 sin( 𝐴𝑜 𝐴) (20)
Finally we use a Newton solver to find the angle of attack that balances the lift coefficient from the aerodynamic surrogate model and the one that allows steady flight equilibrium. In other words we want to solve the following non linear coupled system:

find 𝐴𝑜 𝐴 ∈ [-15; 15] such that 𝐶 𝐿 𝑒𝑞 (𝐶 𝐷 (AoA, Mach, alt)) = 𝐶 𝐿 (AoA, Mach, alt) (21)
The thrust coefficient 𝐶 𝑇 used in Eq. ( 16) is the one of Eq. ( 19) once the Newton solver which balances the lift coefficient has converged. The workflow of the steady flight equilibrium group is described in Fig. 6. Figure 7 shows an optimal trajectory obtained on a climb phase between 2 control points. The phase begins at an altitude of 6 kft and a range of 0 NM and ends at 6.4 kft and 2 NM. The optimal trajectory in terms of energy consumption appears to be linear for this phase. The blue dots representing the DYMOS solution match perfectly with the blue dotted line representing the simulation of the dynamical system evolution when subjected to the controls computed with DYMOS. This match ensures that the optimal control solution is dynamically realistic.

V. Mission test cases

In this section, the overall methodology adopted for the optimization and the mission test cases are described. Four test cases are considered (Mission 1 to Mission 4), each corresponding to a different mission of gradual complexity. Regardless of the test case considered, the optimization problem we want to solve is to maximize the mission end 𝑆𝑂𝐶 with respect to 19 design variables 𝑥 ∈ Ω where the design space 19 and satisfying two structural constraints. The design variables are listed in Table 1 and the constraints in Table 2.

Ω = ( [0, 45] × [-20, 20] × [0, 1.2] × [-5, 5] 3 × [0.001, 0.01] 3 ) 2 ⊂ R
The high-level optimization problem is expressed as follows:

𝑥 * = arg max 𝑥 ∈Ω SOC final (𝑥)
such that wing failure (𝑥) ≤ 0 tail failure (𝑥) ≤ 0 [START_REF] Chauhan | Low-fidelity aerostructural optimization of aircraft wings with a simplified wingbox model using OpenAeroStruct[END_REF] In this paper we want to compare the performance of the Super Efficient Global Optimization (SEGO) algorithm [START_REF] Sasena | Flexibility and efficiency enhancements for constrained global design optimization with kriging approximations[END_REF], a mono fidelity Kriging based Bayesian optimization technique that only uses the high fidelity model with the performance of the Multi-Fidelity Super Efficient Global Optimization (MFSEGO) algorithm [START_REF] Charayron | Multi-fidelity constrained Bayesian optimization, application to drone design[END_REF][START_REF] Meliani | Multi-fidelity efficient global optimization: Methodology and application to airfoil shape design[END_REF], the multi-fidelity co-Kriging based Bayesian optimization method described in Section III.B which takes advantage of both the high and the low fidelity electric drone models. The squared exponential kernel and a constant trend are chosen to build the Kriging (if mono fidelity) or the co-Kriging (if multi-fidelity) models at each iteration of the Bayesian optimizer. The Super Efficient Global Optimization Mixture Of Experts (SEGOMOE) framework [START_REF] Bartoli | An adaptive optimization strategy based on mixture of experts for wing aerodynamic design optimization[END_REF][START_REF] Bartoli | Improvement of efficient global optimization with application to aircraft wing design[END_REF][START_REF] Priem | SEGOMOE: Super Efficient Global Optimization with Mixture of Experts[END_REF][START_REF] Bartoli | Adaptive modeling strategy for constrained global optimization with application to aerodynamic wing design[END_REF] is selected to perform the Bayesian optimizations. SEGOMOE is a Bayesian optimization internal ONERA & ISAE-SUPAERO toolbox. It relies on the Surrogate Modelling Toolbox (SMT) introduced in Section III.A. Note that [START_REF] Vauclin | Développement de modèles réduits multifidélité en vue de l'optimisation de structures aéronautiques[END_REF] implemented MFSEGO in SEGOMOE, [START_REF] Meliani | Multi-fidelity efficient global optimization: Methodology and application to airfoil shape design[END_REF] adapted this implementation to fit the SMT format. It is expected that the use of multi-fidelity will decrease the total optimization cost by reducing the number of calls to the high-fidelity model. One can remark that since the number of design variables is relatively high in our problem it is useful to use a dimension reduction technique in order to save computational cost. To do so we choose to use the partial least squares method (PLS) technique [START_REF] Bouhlel | Improving kriging surrogates of high-dimensional design models by Partial Least Squares dimension reduction[END_REF][START_REF] Ng | A simple explanation of partial least squares[END_REF].

The PLS method finds a linear relationship between input variables and the output variable by projecting input variables onto a new space of lower dimension formed by newly chosen variables called latent variables. PLS has been integrated into the SEGOMOE and can be coupled with SEGOMOE's multi-fidelity kriging method. In our tests, we realized that 4 latent variables were enough to capture most of the information. Adding a fifth latent variable does not improve the precision since almost all of the variance is captured by the 4 first latent variables. On the other hand, removing one of them will strongly reduce the prediction capability without reducing the computational cost that much. Then, we decided to set the dimension of latent space to 4. For the MFSEGO method, the cost ratio between the HF and LF models must be defined. Monte Carlo simulations of the HF and LF models on 100 randomly drawn points according to a uniform distribution are used to estimate the computational costs for all mission test cases. The cost ratios are given by the CPU time ratios. They are summarized in Table 3. The control points defining the four different missions are shown in Fig. 8 where the altitude and Mach values are given as a function of range. In order to compare the performances of the mono and multi-fidelity methods, we performed 10 runs for the single-fidelity Bayesian optimization using SEGO and 10 runs for the multi-fidelity optimization using MFSEGO, respectively. The mono and multi-fidelity runs share the same HF DoE. Table 4 details the sizes of the initial DoE for each mission test case. For each multi-fidelity run, the size of the initial LF DoE is twice the size of the initial HF initial DoE and contains, among other points, the points of the initial HF DoE in order to respect the nested structure of the DoE. The COBYLA [START_REF] Powell | A direct search optimization method that models the objective and constraint functions by linear interpolation[END_REF] optimizer which stands for Constrained Optimization By Linear Approximation is used to determine a reference solution for each test case. COBYLA is an implementation of Powell's nonlinear derivative-free constrained optimization that uses a linear approximation approach. The algorithm is a sequential trust-region algorithm that employs linear approximations to the objective and constraint functions. The obtained optimal designs for each mission test case are shown in Fig. 9. Now that a reference solution is available we can propose a convergence criterion Fig. 9 Radar plot of the optimal drone design for each mission test case.

for our methods. Knowing that a point is considered feasible if the constraints are satisfied at this point, we denote by 𝑦

𝑓 𝑒𝑎𝑠𝑖𝑏𝑙𝑒 𝑏𝑒𝑠𝑡

the best feasible objective value at a given iteration of the Bayesian optimization and ȳ 𝑓 𝑒𝑎𝑠𝑖𝑏𝑙𝑒 𝑏𝑒𝑠𝑡 the respective mean value over the 10 runs. We can define the following convergence criterion method:

| ȳ 𝑓 𝑒𝑎𝑠𝑖𝑏𝑙𝑒 𝑏𝑒𝑠𝑡 -𝑓 HF (ref 𝑠𝑜𝑙 )| ≤ 𝜖 0 | 𝑓 𝐻𝐹 (ref 𝑠𝑜𝑙 )| + 𝜖 1 ( 23
)
where 𝜖 0 is a relative tolerance chosen equal to 0.001 while 𝜖 1 is an absolute tolerance chosen equal to 0.05. The four different missions that the drone must satisfy are detailed in the following sections.

A. Mission 1

Mission 1 is defined by the set of four waypoints detailed in the Table 5. The values in this table are the values imposed on the control points. This first mission consists in a climb phase followed by a cruise phase and a descent phase.

C. Mission 3

Table 7 shows the 8 waypoints that define Mission 3. The flight is performed at low altitude and includes two cruise phases.

D. Mission 4

The 6 waypoints of Mission 4 are detailed in Table 8. The main difference with previous missions is that the cruise phase must be performed at a highest minimal speed. A higher speed is then imposed. This mission consists in a climb phase, an acceleration phase, a cruise phase performed at least at 𝑀𝑎𝑐ℎ = 0.1, a deceleration phase and finally a descent phase.

VI. Results

For the four test cases, the MFSEGO method outperforms the SEGO method. Multi-fidelity appears to improve the budget required for convergence. Figures 10a,10b, 10c and 10d show, respectively for Missions 1, 2, 3, 4, the evolution of the mean and 1𝜎 confidence interval over the 10 runs for the objective value at the best current feasible point. The budget required to achieve the convergence criterion of Eq. ( 23) for both the mono-fidelity and multi-fidelity methods is presented in Table 9 for each mission test case. In this table, "Not reached" means that the convergence criterion was not met with the allowed budget of 100 computational units. One computational unit is the cost of one HF evaluation. For the first 3 missions, Fig. 10a, 10b and 10c show that compared to the SEGO method, the MFSEGO method approaches the optimum faster and eventually reaches the convergence criterion (see Eq. ( 23)) faster too. For Mission 4, at the beginning, both SEGO and MFSEGO methods approach the optimum at the same speed but at the end the convergence criterion (see Eq. ( 23)) is reached with the MFSEGO method and not with SEGO. Finally, for Mission 1, 2 and 3, the use of the bi-fidelity MFSEGO method instead of SEGO reduces the computational budget to reach convergence respectively by a factor 2.24, 2.21, 1.32 respectively. For Mission 4, the mono-fidelity method does not even converge for each run. Only 6 runs out of 10 reached the optimum with the allowed budget in the case of the mono-fidelity method, while all runs reached the optimum in the case of the multi-fidelity method. The multi-fidelity method is therefore much more robust.

VII. Conclusion

In this paper, we addressed the problem of designing an electric fixed-wing drone for long range surveillance missions. We recall that the methodology applied to solve the optimization problem relies on a multi-fidelity Bayesian optimization. We developed two drone models of different fidelities. The main contribution here is the addition of an electric propulsion module and of a mission module to these drone models. For the low fidelity drone model, a mission component based on the electrical Breguet range equation has been developed while in the high fidelity model, the mission discipline has been tackled using DYMOS, performing optimal control on the mission trajectory. The use of metamodels to approximate the lift and drag coefficients provides sufficient accuracy while greatly reducing the cost of calling Dymos since we no longer have to call OpenAeroStruct many times when solving the optimal control problem. In the end we used both of these models to perform mono and multi-fidelity Bayesian optimizations on four different missions that lead to different drone designs. Finally, to solve the optimization problem, the MFSEGO multi-fidelity algorithm overtakes SEGO mono-fidelity algorithm in terms of needed convergence budget and robustness.

B. Classic fuel Breguet range equation

The classic fuel Breguet range equation is given by:

𝑅 = - ℎ 𝑔 𝐿 𝐷 𝜂 overall ln 𝑊 final 𝑊 init ( 24
)
where ℎ is the fuel energy per unit of mass, 𝑔 is the gravity acceleration, 𝐿 is the lift force, 𝐷 the drag force, 𝜂 overall is the overall efficiency (ie: propulsive power over fuel power), 𝑊 final is the aircraft final weight and 𝑊 init is the initial aircraft weight.

C. Electric Breguet range equation demonstration

First note that the range 𝑅 can be expressed with the speed and the time flying at that speed: 𝑅 = 𝑉 × 𝑡. Under ideal conditions, the time needed to empty a battery is given by: 𝑡 = 𝑚 batt 𝑒 𝑏 𝑃 batt [START_REF] Priem | SEGOMOE: Super Efficient Global Optimization with Mixture of Experts[END_REF] where 𝑒 𝑏 is the specific energy density of the battery (unit: J.kg -1 ), 𝑚 batt is the battery mass and 𝑃 batt is the battery power. Lets define the overall propulsion system efficiency 𝜂 𝑡𝑜𝑡 = 𝑃 prop 𝑃 batt with 𝑃 prop is the required propulsive power to

Note that the 𝐿 𝐷 ratio can be replaced by the ratio of lift and drag coefficient 𝐶 𝐿 𝐶 𝐷 . Indeed, 𝐿 = 𝑞𝑆𝐶 𝐿 and 𝐷 = 𝑞𝑆𝐶 𝐷 with 𝑆 being the lifting surface area and 𝑞 being the dynamic pressure.

  𝐶 𝐷 = Drag coefficient 𝐶 𝐿 = Lift coefficient 𝐷 = Drag force (N) 𝑒 𝑏 = Battery specific energy (Wh.kg -1 ) 𝐸 𝑆𝐶 = Electronic speed controller 𝜂 𝐸𝑆𝐶 = ESC efficiency 𝜂 𝑚 = motor efficiency 𝜂 𝑝 = propeller efficiency 𝜂 𝑡𝑜𝑡 = total efficiency 𝑔 = Gravity acceleration (m.s -2 ) 𝛾 = Flight path angle 𝐿 = Lift force (N) 𝑚 = Total drone mass (kg) 𝑚 batt = Battery mass (kg) 𝑀 𝐷 𝐴 = Multidisciplinary Analysis 𝑀 𝐷𝑂 = Multidisciplinary optimization 𝑂𝐷𝐸 = Ordinary differential equation 𝑃 batt = Battery output power 𝑃 𝐸𝑆𝐶 = ESC output power 𝑃 𝑚 = Motor output power 𝑃 prop = Propeller output power 𝑄 batt = Capacity of the battery (Coulomb) 𝑆𝑂𝐶 = State of charge (%) 𝑇 = Thrust force (N) 𝑇 𝐴𝑆 = True airspeed (s) 𝑈 𝐴𝑉 = Unmanned aerial vehicle 𝑈 batt = Total electric tension of the battery (V) 𝑈 cell = Electric tension of a single cell of the Lithium-Polymer (LiPo) battery (V) 𝑊 = Weight force (N)

Fig. 2

 2 Fig. 2 XDSM of the RootChordTail group.

Fig. 5

 5 Fig. 5 XDSM of the propulsion group.

Fig. 6

 6 Fig. 6 XDSM of the steady flight equilibrium group which is a sub analysis of Fig. (6).

Fig. 7

 7 Fig. 7 Illustration of an optimal trajectory obtained on a climb phase.

  control points along wing span [0.001, 0.01] 3 m 3 structural tube thickness control points along tail span [0.001, 0.01] 3 m 3 twist control points along wing span [1, 1.5]𝑥 [0.5, 1]𝑥 [0, 0.5] degree 3 twist control points along tail span [1, 1.5]𝑥 [0.5, 1]𝑥 [0, 0.5] degree 3 Total number of design variables 19

Fig. 8

 8 Fig. 8 Range, altitude (red dot) and Mach values (blue cross) at the control points of the four missions.

Fig. 10

 10 Fig. 10 Mean and 1-sigma confidence interval over 10 runs of the objective value at the best HF feasible point against budget.

Fig. 11 MFSEGO

 11 Fig. 11 MFSEGO methodology diagram.

Table 1 Definition of the 19 design variables.

 1

	variables	design space	unit	quantity
	wing span	[5, 6]	m	1
	root chord wing	[0.9, 1.2]	m	1
	wing taper ratio	[0.6, 1]	no unit	1
	wing sweep angle

Table 2 Definition of the constraints.

 2

	name	type	bound	dimension tolerance
	wing failure	≤	0	1	1𝑒 -7
	tail failure	≤	0	1	1𝑒 -7
	Total number of constraints	2

Table 3 Cost ratio for each mission test case.

 3

		Mission 1 Mission 2 Mission 3 Mission 4
	cost ratio	140	124	168	164

Table 4 Initial DoE size for each mission test case.

 4

		Mission 1 Mission 2 Mission 3 Mission 4
	HF	10	10	5	5
	LF	20	20	10	10

Table 5 Range, altitude, Mach and climb rate values at the Mission 1 waypoints.

 5 Mission 2 consists of by six waypoints (see Table6). The specificity of this mission is that it carries out a level climb. Moreover, the flight altitude is much lower.

	Range	alt	Mach Climb rate
	control point 0 0 NM 6.0 kft
	control point 1 2 NM 6.4 kft 0.065	0 ft min
	control point 2 7 NM 6.4 kft 0.065	0 ft min
	control point 3 9 NM 6.0 kft
	B. Mission 2

Table 6 Range, altitude, Mach and climb rate values at the Mission 2 waypoints.

 6

Table 7 Range, altitude, Mach and climb rate values at the Mission 3 waypoints.

 7

	Range	alt	Mach Climb rate
	control point 0 0 NM 1.0 kft
	control point 1 2 NM 1.4 kft 0.065	0 ft min
	control point 2 6 NM 1.4 kft 0.065	0 ft min
	control point 3 7 NM 1.2 kft 0.065	0 ft min
	control point 4 9 NM 1.2 kft 0.065	0 ft min
	control point 5 11 NM 1.6 kft 0.065	0 ft min
	control point 6 13 NM 1.6 kft 0.065	0 ft min
	control point 7 16 NM 1.0 kft

Table 8 Range, altitude, Mach and climb rate values at the Mission 4 waypoints.

 8

	Range	alt	Mach Climb rate
	control point 0 0 NM 6.0 kft
	control point 1 2 NM 6.4 kft 0.065	0 ft min
	control point 2 3 NM 6.4 kft 0.12	0 ft min
	control point 3 10 NM 6.4 kft 0.12	0 ft min
	control point 4 11 NM 6.4 kft 0.065	0 ft min
	control point 5 13 NM 6.0 kft

Table 9 Needed budget to reach convergence criterion for each method and each mission test case.

 9

		Mission 1 Mission 2 Mission 3	Mission 4
	Mono-fidelity (1F)	30	39	29	Not reached
	Multi-fidelity (2F)	13.39	17.65	21.90	30.61

Acknowledgments

The PhD is funded by the defense innovation agency (AID) and by the french Directorate General of Armaments (DGA) as part of the CONCORDE project. This work is part of the activities of ONERA -ISAE -ENAC joint research group.

Appendix

In this appendix, Section VII.A presents the MFSEGO methodology in a diagrammatic form. Then, the classical Breguet fuel autonomy equation is recalled in Section VII.B and Section VII.C details how to find the expression of the electric Breguet range equation.