Yannick Schwartz
email: yannick.schwartz@inria.fr

Gaël Varoquaux

Christophe Pallier

Philippe Pinel

Jean-Baptiste Poline

Bertrand Thirion

Keywords: Meta-analysis, fMRI, multiple comparison, machine learning

Typical cohorts in brain imaging studies are not large enough for systematic testing of all the information contained in the images. To build testable working hypotheses, investigators thus rely on analysis of previous work, sometimes formalized in a so-called meta-analysis. In brain imaging, this approach underlies the specification of regions of interest (ROIs) that are usually selected on the basis of the coordinates of previously detected effects. In this paper, we propose to use a database of images, rather than coordinates, and frame the problem as transfer learning: learning a discriminant model on a reference task to apply it to a different but related new task. To facilitate statistical analysis of small cohorts, we use a sparse discriminant model that selects predictive voxels on the reference task and thus provides a principled procedure to define ROIs. The benefits of our approach are twofold. First it uses the reference database for prediction, i.e. to provide potential biomarkers in a clinical setting. Second it increases statistical power on the new task. We demonstrate on a set of 18 pairs of functional MRI experimental conditions that our approach gives good prediction. In addition, on a specific transfer situation involving different scanners at different locations, we show that voxel selection based on transfer learning leads to higher detection power on small cohorts.

Improving accuracy and power with transfer learning using a meta-analytic database

Multi-subject or multi-condition experiments are the workhorse of bio-medical imaging research, whether it be drug development or basic research. Imaging provides a wealth of information on the biomedical problem at hand. However the typical sample size is too small to fully exploit this information. For this reason, investigators often turn to previous studies in order to formulate hypotheses and restrict the search space, i.e. select a subset of anatomical or functional structures of interest to the current study. A typical case is that of early-stage clinical trials, for which the group size is very small, but that are most often based on previous results concerning the pathology under study. However, understanding the literature is increasingly difficult and requires a systematic approach, that takes the form of a meta-analysis, pooling results from multiple experiments that address a set of related research hypotheses [START_REF] Sutton | Methods for metaanalysis in medical research[END_REF]. In particular, brain imaging studies heavily rely on such meta-analyses [START_REF] Wager | Evaluating the consistency and specificity of neuroimaging data using meta-analysis[END_REF], as the brain is still an ill-understood and complex organ. In functional Magnetic Resonance Imaging (fMRI) studies, typical group sizes range from 10 to 20 subjects, which is not always enough to warrant the reliability of brain-wide analysis [START_REF] Thirion | Analysis of a large fmri cohort: Statistical and methodological issues for group analyses[END_REF]. More importantly, the time that can be spent in the scanner by subjects is limited, and not all interesting experimental conditions will be acquired.For this reason, it is common practice to reduce the study to a set of regions of interest (ROIs) extracted from the literature. Investigators define these ROIs by extracting locations of peak activations from the literature [START_REF] Yarkoni | Largescale automated synthesis of human functional neuroimaging data[END_REF], or from coordinate databases such as BrainMap [START_REF] Laird | Brainmap: the social evolution of a human brain mapping database[END_REF]. While most of these meta-analyses are conducted on activation coordinates, the increase of data sharing opens the door to meta-analysis on full brain images which results in higher statistical power [START_REF] Salimi-Khorshidi | Metaanalysis of neuroimaging data: a comparison of image-based and coordinate-based pooling of studies[END_REF]. Previous statistical and modeling work on meta-analysis for fMRI has focused on better modeling of the reference database [START_REF] Wager | Evaluating the consistency and specificity of neuroimaging data using meta-analysis[END_REF].

In this work, we are interested in the generalization power of meta-analyses on new data. We introduce a new meta-analysis method using a reference database of images to guide statistical analysis of a new dataset. In particular we rely on predictive models, useful to learn biomarkers, and use them to select relevant voxels in order to increase the statistical power of a new study.

Methods

Problem setting We start from a reference database made of l experiments, each comprising n l contrasts possibly acquired over multiple subjects. We denote the brain images by X l ∈ R n l ,p with associated experimental condition y l ∈ R n l . Given a new experiment, denoted target, (X , y ) ∈ (R n ,p , R n ), we are interested in two problems: i) (biomarkers) can we predict y from X ? ii) (inference) can we test hypotheses on the links between y and X , for instance in a linear model? These are ill-posed problems from the statistics standpoint, as n p. The root of the problem is the dimensionality of the data: medical images are composed of many voxels, typically p ≈ 50 000 with fMRI. This large number of descriptors limits statistical inference power due to multiple testing; a problem that appears in predictive approaches as the curse of dimensionality. Here, we use our reference database to better condition this statistical problem.

Transfer learning The gist of our approach is to learn on some experiments of our database (X l , y l ) discriminative models that contain predictive information for the target experiment (X , y ). In machine learning, this problem is known as transfer learning [START_REF] Pan | A survey on transfer learning[END_REF]. The underlying assumption of transfer learning is the same as that for meta-analysis: the reference database should contain some common information with the target experiment. Here we use a simple form of transfer learning: we train a linear classifier on an experiment in the database that is similar from the neuroscientific point of view to the new data, and use it to predict the labels of the new data.

Selecting predictive features We use a sparse linear classifier, specifically an 1penalized logistic regression. The motivation behind this choice of classifier is that it produces a sparse set of weights that can be used to select relevant voxels. In particular, under certain conditions, the classifier can recover with high probability the complete set of k features in X that are predictive of y for a sample size of n min = O k log p [START_REF] Bach | Self-concordant analysis for logistic regression[END_REF]. The logarithmic dependence in p is an appealing property in view of the dimensionality of medical imaging datasets.

In practical situations, it can be hard to control the errors on this feature selection, in particular as it depends on the choice of the amount of 1 penalty. For this reason, Meinshausen and Bühlmann [START_REF] Meinshausen | Stability selection[END_REF] introduce randomized variants of sparse estimators, that can be seen as sampling the posterior probability of selection and keeping only features that are selected frequently. In particular, they establish non-asymptotic recovery results for the randomized lasso, which consists in applying the Lasso on random subsamples of the data and rescaling of the regressors. Here, we adapt this strategy to classification as the logistic regression is locally equivalent to a weighted least square and recovery results can carry from square-loss regression to logistic regression [START_REF] Bach | Self-concordant analysis for logistic regression[END_REF].

We want to use transfer learning to perform screening of the voxels, i.e. eliminate many voxels that are not related to our target experiment. For this purpose, we need a low probability of rejecting relevant variables. Each iteration of the sparse logistic regression in the randomized logistic can select reliably only k max ≈ n/ log p variables. In the worst case situation, we have k heavilycorrelated variables and one of them is selected at random by the sparse logistic regression at each iteration. For each of these variables, the probability of selecting it less than s times during m iterations of the randomized logistic is given by the cumulative distribution function of a binomial with per trial success ratio 1/k. If s ≤ m/k, by Hoeffding's inequality, this probability goes to zero in o exp m . In other words, if we impose a threshold τ = s/m on the selection frequency, we can recover a group of k correlated variables for τ ≤ 1/k.

Brain parcellations

Although randomization relaxes the conditions on recovery, a remaining necessary condition is that the regressors of interest, i.e. the values x i across the subjects of the k predictive voxels, must be weakly correlated 4 . Because of the large amount of smoothness present in medical images, in particular in group-level fMRI contrasts, these conditions cannot be satisfied. Indeed, values taken by a voxel are very similar to values taken by its neighbors. In addition, the numbers of subjects used in fMRI are often below the sample size required for good recovery. For these reasons we resort to feature agglomeration: using hierarchical clustering to merge neighboring voxels carrying similar information into parcels [START_REF] Michel | A supervised clustering approach for fmri-based inference of brain states[END_REF]. This strategy brings the double benefit of reducing the problem size, and thus the required sample size, and mitigating local correlation, at the expense of spatial resolution.

Experiments and Results

FRMI datasets

We use 3 studies for this meta-analysis. The first study (E1 ) [START_REF] Pinel | Fast reproducible identification and large-scale databasing of individual functional cognitive networks[END_REF] is composed of 322 subjects and was designed to assess the inter-subject variability in some language, calculation, and sensorimotor tasks. The second study (E2 ) is similar to the first one in terms of stimuli, but its data was acquired on 35 pairs of twinsubjects. The last study (E3 ) [START_REF] Pallier | Cortical representation of the constituent structure of sentences[END_REF] characterizes brain regions in charge of the syntactic and the semantic processing for the language. It was performed with 40 subjects, 20 of which were stimulated by pseudowords (jabberwocky stimuli) instead of actual meaningful sentences. All the studies were pre-processed and analyzed with the standard fMRI analysis software SPM5. The data used for this work are the statistical images resulting from the intra-subject analyses across the 3 studies. E1 has 34 contrasts images available, E2 56, and E3 28. The raw images were not always acquired on the same scanner. E1 has data from a 3T SIEMENS Trio, and a 3T Brucker scanner; E2 data were acquired on a 1.5T GE Signa; and E3 images come from the same 3T SIEMENS Trio.

Experimental results for prediction

Here we are interested in the prediction problem: using transfer learning to discriminate a pair of constrasts with an estimator trained on two other contrasts.

We used 4 different approaches to learn the discriminative models. The first approach relies on the activation likelihood estimate (ALE) method [START_REF] Laird | Ale meta-analysis: controlling the false discovery rate and performing statistical contrasts[END_REF], as this is a commonly published method for coordinate-based meta-analyses. We extract the activation positions from the contrasts maps, and then apply a Gaussian kernel. We use the preferred FWHM of 10mm [START_REF] Turkeltaub | Meta-analysis of the functional neuroanatomy of single-word reading: method and validation[END_REF]. The other approaches directly use the contrast images. We name raw contrasts the method based on the contrasts voxels values; contrast-specific parcels the method that uses parcels from the training set: and meta-analystic parcels the method that learns parcels from the full database. We evaluate on our base of contrasts the ability to do transfer learning, i.e to learn decision rules that carry over from one situation to another. Since we must make the assumption that the reference contrasts hold common information with the contrasts of interest, we do not try out all the possible combinations, but rather manually select pairs of contrasts from a single experiment that form a meaningful classification task (e.g., computation versus reading, or Korean language versus French language). Out of all the possible combinations, we select 35 pairs of classification task, and subsequently combine them into 18 transfer pairs, on which it is reasonable to think that the transfer could occur (e.g., computation and reading in visual instructions, transfer on computation and reading in auditory instructions). We first train a linear classifier within 6-fold cross validation test on a first set of pairs, setting the penalization amount by nested cross-validation, we call this step inline learning. We then re-use the discriminative model on a different pair of contrasts to perform the transfer learning. The 3 studies containing language related tasks, we can transfer between pairs within an experiment, and across experiments. Among the 18 selected transfer pairs, we find that 9 can give rise to such a transfer. Since a transfer is directed, we perform it both ways, which yields once again 18 transfer pairs to test upon. The associated prediction scores from the different methods are reported in Table 1. The general observation is that ALE yields a poorer prediction performance than any other method. This is true both for the transfer and inline predictions. We also find that brain parcellations scores similar to the raw contrasts images, and closer to the inline predictions. We find that while the contrast-specific parcels and meta-analytic parcels methods do not return the same parcels, they produce very close results. We can thus use the full database to learn a single reference parcellation to perform meta-analysis.

Experimental results for inference

Here we are interested in the inference problem: using transfer learning to help hypothesis testing on a target dataset. In the following, we only consider a specific transfer, namely the last line in Table 1: we learn a model discriminating French native speakers reading French or Korean, and apply it on another experiment in which French subjects had to read French or jabberwocky. This transfer is interesting as it involves two different experiments acquired on different scanners, and cognitive paradigms that share a similar expression, incomprehension of visual language stimuli. As can be seen in Table 1, the prediction scores of transfer learning as well as inline learning on this pair are acceptable although not excellent: French language and jabberwocky are difficult to separate.

Figure 2 gives the stability scores of the randomized logistic discriminating reading Korean from reading French for the different set of features -activation peaks, raw contrasts, parcels learned on the training contrasts or on the full database. We can see that while learning at the voxel level or at the parcel level gives similar prediction performance (Table 1), the stability score maps are very different. At the voxel-level, with 70 subjects (p = 40 000, n = 70) the recovery is limited to approximately 7 voxels without randomization: the recovery conditions are violated. As a result, the randomized logistic selects only the most predictive voxels. On the parcels, contrast-specific or meta-analytic (i.e., learned on the full database), the selection frequency highlights regions of the brain that are known to be relevant for language comprehension, including the left anterior superior temporal sulcus and the part of the temporal parietal junction (Wernicke's area).

We threshold the stability selection scores of the first experiment (Korean vs French) to select candidate voxels for the target experiment (jabberwocky vs French). As we want to perform a rough screening and would rather err on the side of false detections than false rejections, we take a very low threshold τ = .01. Following our analysis above, the size of the largest group of correlated features that we can detect with such a threshold is on the order of 1/τ ≈ 100. With 2000 parcels, this number corresponds to 5% of the brain, i.e. 8 000 voxels, and we can safely consider that no fMRI contrasts is composed of groups of heavily correlated features larger than this fraction.

On the target experiment, we perform a standard group-level analysis with the voxels selected, testing for the difference between the two conditions, jabberwocky or French reading. We report results with p-values corrected for multiple comparisons at a given family-wise error rate (FWER) using Bonferroni correction, and for a given false discovery rate (FDR) using the Benjamini-Hochberg procedure. On table 2, we compare the number of detections and the detection rate, i.e. the fraction of voxels detected as significantly different, for a full brain analysis and for an analysis limited to the voxel selection. We compare our voxel selection method to a one-way ANOVA, and find that transfer learning outperforms the ANOVA for all the cohort sizes. Figure 3 shows the Q-Q plots on which the Benjamini-Hochberg procedure is applied. We find that voxel selection by transfer learning improves both the absolute number of detections and the detection rate for FWER and FDR correction.

Conclusion

In this paper, we propose to improve the conditioning and power of statistical analyses in imaging studies, using a large meta-analytic database.

In a transfer learning scheme, we train on the database sparse discriminative models that are suited to the target experiment. Not only can the predictive power of these models can be useful to establish biomarkers, but also they perform feature selection that can increase the statistical power of a standard group analysis on new experiments, provided enough predictive features (voxels) can be recovered. Using brain parcellations, the discriminative model acts to screen parcels unlikely to be relevant in the target experiment, thus defining automatically ROIs.

Using a set of 3 fMRI studies related to language, we confirm experimentally that our transfer learning scheme is able to: i) perform accurate predictions on experiments acquired on a different scanner and with varying paradigm, ii) outperform the standard meta-analysis procedures based activation peaks, iii) increase the statistical power in the target experiment by using the ROIs defined by the discriminative model.

In this work we manually select the contrast pairs since it is delicate to interpret a transfer learning score without good knowledge of the cognitive or clinical conditions under study. Future work will study automatic contrast pairs selection, e.g. by mining the descriptions of the experiments [START_REF] Yarkoni | Largescale automated synthesis of human functional neuroimaging data[END_REF], to address the problem of synthesizing the ever-growing literature and data in medical research.

Fig. 1 .Fig. 2 .Fig. 3 .

 123 Fig. 1. Prediction performance relative to the best performing approach: inline prediction with raw contrasts images: the p-values indicate whether the associated methods are significantly poorer than the best performing method. Activation Peaks Contrast-specific Parcels Meta-analytic Parcels Raw contrasts 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0.0 0.1 Relative prediction performance

Table 1 .

 1 Prediction scores for inline and transfer learning. trans.= transfer; in.= inline; comp.= computation, sent.= sentences (reading), jabb.= jabberwocky; S= sentence with one word constituents, L= one constituent long sentence.

	Names	Peaks Contrasts Parcels Meta parcels
		trans. in. trans. in. trans. in. trans. in.
	E1, comp./sent. → E2, comp./sent.	0.75 0.85 0.88 0.97 0.83 0.96 0.83 0.96
	E2, comp./sent. → E1, comp./sent.	0.66 0.83 0.88 0.96 0.85 0.95 0.85 0.96
	E3, jabb./French (L) → E3, jabb./French (S) 0.46 0.48 0.65 0.67 0.62 0.60 0.67 0.62
	E3, jabb./French (S) → E3, jabb./French (L) 0.52 0.71 0.67 0.85 0.71 0.85 0.65 0.79
	E3, jabb./French (L) → E2, Korean/French 0.65 0.46 0.73 0.79 0.65 0.81 0.76 0.85
	E2, Korean/French → E3, jabb./French (L) 0.73 0.81 0.79 0.85 0.75 0.81 0.75 0.75

Specifically, the condition for recovery with randomized lasso it is a lower bound on the conditioning of the sparse eigenvalues of the design matrix [9, theorem 2] and for sparse logistic regression the corresponding condition is a lower bound on the eigenvalues of the regressors of interest's covariance matrix[8, theorem 4].

Acknowledgements This work was supported by the ANR grants BrainPedia ANR-10-JCJC 1408-01 and IRMGroup ANR-10-BLAN-0126-02.

FWER corrected FDR corrected n All voxels

Selection ANOVA All voxels Selection ANOVA 10 0 (0%) 0 (0%) 0 (0%) 0 (0%) 0 (0%) 0 (0%) 20 0 (0%) 3 (0.02%) 0 (0%) 0 (0%) 4 (0.027%) 0 (0%) 40 5 (0.0084%) 33 (0.22%) 2 (0.0014%) 143 (0.97%) 1339 (9%) 201 (1.4%) Table 2. Number of detections at p < 0.05 for difference cohort size, for transfer learning and ANOVA. The percentage of detection is indicated in parenthesis.


Coquelin L Fischer

Le Brusquet

Aerosol size distribution measurement using a SMPS: scanning mode and uncertainty analysis

Keywords: Aerosol Size distribution measurement, Electrical classification, uncertainty analysis, Data inversion

A model to simulate SMPS (Scanning Mobility Particle Sizer) measurement and the associated uncertainty analysis when axial DMA (Differential Mobility Analyser) classifier operates under scanning mode conditions is described. Starting from simulated SMPS raw data, a fast estimation of aerosol size distribution measurement using regularization technique is performed. Then, global sensitivity analysis is used to discriminate significant parameters of the system and, as a preliminary result, a 95% confidence region is obtained by Monte Carlo simulations on an atmospheric aerosol size distribution.

I. Introduction

A MONG devices using electrical mobility methods to measure aerosol particle size distribution, the SMPS (Scanning Mobility Particle Sizer) is the most widely used to characterize particles smaller than 0.1µm. Since it can be applied in fields like: indoors air quality measurement, vehicle exhaust pipes, atmospheric studies, toxicology testing, etc., it has become of interest to assess the uncertainty associated with its outputs. So far, although some works have been carried on the adequacy of particle size measurements made by this device, such a dimensional feature uncertainty cannot be provided. First, a full SMPS model is presented when DMA operates under scanning mode. Then, an alternative to common data inversion performed in aerosol science that allies fast reconstruction of the aerosol size distribution and accuracy is described. Fast reconstruction of the solution is a requirement here due to the fact that when dealing with uncertainty propagation in a numerical code, a large number of calls is necessary. Last section pays a particular attention to uncertainty analysis and Monte-Carlo simulations lead to a 95% confidence region on an atmospheric aerosol size distribution. Figure 1 illustrates how to perform one simulation of the numerical code: starting from a known aerosol size distribution n(d p ), d p ∈ [d min , d max ] being the interval of scanned diameters, the code returns an estimate n(d p ).

II. Scanning mode theory

SMPS is composed of a Differential Mobility Analyser (DMA) and a Condensation Particle Counter (CPC). First one selects particles mobility diameters while the other is an optical counter used for the detection phase. First, aerosol passes through an inertial impactor to avoid largest particles to enter the DMA column, then the aerosol enters the charge neutraliser to be conditioned. Indeed, aerosol particles are usually charged. TSI neutraliser is used here to create a bipolar charge distribution by exposing aerosol to a radioactive Krypton source, so, particles that carry several charges lose their charge excess. Once aerosol is well conditioned particles are selected using electrical classification inside DMA column: an electric field is created, and the airborne particles drift in the DMA according to their electrical mobility. Flow rates system inside the DMA column and also geometrical parameters are described in the figure 2. L DM A , r 1 and r 2 represent the length of the column, radii of the inner and outer electrodes, respectively. Q sh is the clean sheath air flow, Q a being the aerosol flow, Q cpc the mono disperse flow that goes to the aerosol detector and finally Q exc the excess flow.

Once selected, particles reach the CPC to be counted. Scanning mode is based on Wang and Flagan research [START_REF] Wang | Scanning electrical mobility spectrometer[END_REF]: voltage at the center rod of the DMA is ramped exponentially to ensure particles trajectories to remain unchanged during the all scanning time. Let t s the scanning time, V min and V max respectively the minimum and maximum value of the voltage defined by the user, Collins and al. [START_REF] Collins | The scanning dma transfer function[END_REF] give the following expression for the tension:

V (t) = V min exp(t/τ ) , with τ = t s / ln (V min /V max ) .
CPC detects particles carried by flow rate Q cpc during time ∆ t . Usually, CPC results are expressed in terms of concentration. Choice is made here to work with CPC raw data. Let d p be the particle mobility diameter, Z d the corresponding electrical mobility, p the number of charges carried by each particle, ∆φ(V ) a function of the geometry and the supplied voltage of the DMA and n(d p ) the aerosol size distribution, then basic equation for electrical mobility measurement leads to the CPC raw data response:

r(t + ∆ t ) = +∞ 0 k d p , Z d (p, d p ), ∆φ(V (t)) n(d p )d dp (1)
where r(t + ∆ t ) is the number of particles counted by the CPC at time t + ∆ t for particles selected in the DMA column set at voltage V (t) that is the tension of selection defined as the mean tension for a passage inside the column of the classifier: Indeed, particles reach the CPC only after transport through the DMA and transport through the tubing which connects the CPC with the mono disperse outlet port of the DMA. Transport times are respectively called t f and t d . First is assumed to be equal to the mean residence time of the aerosol inside the DMA column and second depends on the internal CPC delay time and flow rates system (default value is for 10 inches connection to CPC in our case), so

V (t) = 1 t f t-t d t-t d -t f V (ξ)
t f = π(r 2 2 -r 2 1 )L DM A Q sh + Q a ,
DMA and CPC are respectively characterized by the transfer function Ω and size detection efficiency n CP C when charge conditioner is modelled with function Ψ. For a given time t, the kernel function k is given by the following expression:

k d p , Z d (p, d p ), ∆φV (t) = Q a ∆ t +∞ p=1 Ψ(p, d p )n CP C (d p )Ω Z d (p, d p ), ∆φV (t) .
Charge conditioning model is extracted from both Wiedensohler [START_REF] Wiedensohler | An approximation of the bipolar charge distribution for particles in the submicron size range[END_REF] and Gunn [START_REF] Gunn | Diffusion charging of atmospheric droplets by ions, and the resulting combination[END_REF] while n CP C is obtained by fitting experimental data made by Quant [START_REF] Quant | Performance of condensation particle counters with three continuous-flow designs[END_REF] by using cubic interpolation.

Next, modelling DMA transfer function requires the electrical mobility to be known for each mobility diameter d p . Let e the charge of the electron, µ g the dynamic gas viscosity, C c the slip correction factor, electrical mobility is then linked to a specific mobility diameter d p via:

Z d (p, d p ) = (peC c (d p ))/(3πµ g d p ), (2)
where C c is defined as

C c (d p ) = 1 + K n (d p ) [α + β exp(-γ/K n (d p ))] , K n (d p ) = (2λ m )/d p ,
where K n is the Knudsen number and λ m is the mean free path of a particle. Triplet (α, β, γ) is taken from experiments made by Kim and al. [START_REF] Kim | Slip correction measurements of certified psl nanoparticles using a nanometer differential mobility analyser (nano-dma) for knudsen number from 0.5 to 83[END_REF].

An other expression of the electrical mobility, denoted by Z v is:

Z v = ν r /E, (3)
ν r being the drift velocity and E the electric field strength.

Starting from equation ( 3), solving electromagnetic wave equation while assuming: quasi stationary system, uniformity of the velocity profiles inside the DMA column, and laminar flow regime instantaneously reached, one can get the electrical mobility Z v for a set voltage V at the center rod voltage of the DMA:

Z v (V ) = Q 4π∆φ(V ) , with ∆φ(V ) = L DM A V ln(r 2 /r 1 ) . (4)
From both expressions (2) and ( 4), knowing the voltage at which the particle has been selected, one can get its corresponding mobility diameter. Neglecting diffusional effects inside the column of the classifier, Whitby [START_REF] Knutson | Accurate measurement of aerosol electric mobility moments[END_REF] showed the following expression for the DMA transfer function:

Ω(Z d , ∆φ(V (t))) = 1 Q a max(0, [min(Q a , Q cpc , [(Q a + Q cpc )/2 -|2π Z d ∆φ(V (t)) + (Q sh + Q exc )/2|])]).
Until there, a deterministic model that takes into account the physical considerations of the system has been built. The physic is here represented through analytical functions, so the need of having a code that runs quickly is fulfilled.

III. Data Inversion

Reconstruction of an unknown aerosol size distribution is a general mathematical problem: starting from CPC measurements, problem consists in the estimation of original concentration n. Through years several techniques have been found to reach this purpose [START_REF] Seinfeld | Inversion aerosol size distribution data[END_REF], [START_REF] Hagen | Linear inversion method to obtain aerosol size distributions from measurements with a differential mobility analyser[END_REF], [START_REF] Twomey | Comparison of constrained linear inversion and an iterative non-linear algorithm applied to the indirect estimation of particle size distributions[END_REF], etc.. Nevertheless, industrial approaches have to fulfill other constraints such as a fast reconstruction of the signal. Main problem here is to reach enough accuracy in a short time.

A. Common approach

Since Spectrometers operate at high resolution, at each time of the scanning process, the support of Ω is narrow and it is common to consider that only a particle of diameter d * p (t) is extracted from the classifier. When assuming that particles only carry one charge, this diameter is computed by solving the following equation:

Z d (p = 1, d * p (t)) = Z v (V (t)
), and if n denotes the estimate of the aerosol size distribution, its expression takes the following form:

n(d * p (t)) = r(t + ∆ t ) (∂Z d (1, d p )/∂d p ) dp=d * p (t) Q a ∆ t Ψ(1, d * p (t))n CP C (d * p (t)) Z + v Z - v Ωd Z d
, where

Z ± v = Z v (V (t)) ± ∆Z v (V (t)) = 1 ± Q a + Q cpc Q sh + Q exc Z(V (t)).
Multiple

B. Alternative approach

Estimation of the aerosol size distribution can be performed regardless of the DMA resolution being used and without assumptions concerning the number of charges carried by each particle. The method proposed here is based on the discretization of the mobility diameter space. It is somehow comparable to what has been developed by Alofs and Hagen [START_REF] Hagen | Linear inversion method to obtain aerosol size distributions from measurements with a differential mobility analyser[END_REF] when the DMA operates under stepping mode.

First step is to discretize the mobility diameter space using log scale. Let (d p,j ) j=1,..,J the points of reconstruction of the solution, because the particles diameters typically range over several orders of magnitude, a convenient choice for integrating (1) is the variable transformation f (d p ) = ln(d p ). Second step is to foresee the behaviour of the size distribution over a narrow diameter range. Linear interpolation is chosen here because changes in the solution are slow. Moreover, the more points is used to interpolate, the less the assumption will be strong. When using notations f = f (d p ) and f j = f (d p,j ), equation ( 1) is further converted to:

r(t + ∆ t ) = f J f1 d p k f, Z d (p, f ), ∆φV (t) n(f )d f ,
and it can be written in a matrix form as R = KN , where both solution vector N and vector of the observations R are defined as

R i = r(t i ) , t i = t 0 + i∆ t , N j = n(f j )
, j = 1, ..., J, while element (i, j) of matrix K corresponds to the following sum of integrals:

fj fj-1 d p k(f )L j-1 (f )d f + fj+1 fj d p k(f )(1 -L j (f ))d f ,
where L j is the linear interpolant. Good point of this method is that multiple charges correction won't have to be performed any longer because it is corrected by construction.

Knowing R and K expressions, problem now is to retrieve N . Nevertheless, K is badly conditioned. Regularisation technique described by Tikkhonov [START_REF] Tikhonov | Solutions of ill-posed problems[END_REF] is used to overpass this problem. So, if typical least-square method leads to the minimization of the quadratic error KN -R 2 , Tikkhonov suggests to minimize KN -R 2 +λ N 2 . A prior on the solution to rebuild can be carried by the second term multiplied by regularisation parameter λ. Another prior is the positivity of the searched vector N as it represents a concentration. As decided by Seinfeld [START_REF] Seinfeld | Inversion aerosol size distribution data[END_REF], non smooth solutions are penalised here and the minimization becomes min

N ≥0 KN -R 2 + λ DN 2 ,
where D is the finite difference matrix chosen to approximate the integral of the second derivative of N :

D =       -1 2 -1 0 • • • 0 0 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 0 0 • • • 0 -1 2 -1       .
Phillips-Twomey method is used here: a regularization parameter λ is chosen to calculate the corresponding regularization matrix K, K = K t K + λD t D and then minimization of KN -R 2 using Non-Negative-Least-Squares developed by Lawson [START_REF] Lawson | Solving least squares problems[END_REF] is performed. The regularization hyperparameter is chosen with a graphical technique using the L-curve criterion [START_REF] Hansen | Analysis of discrete ill-posed problems by means of the l-curve[END_REF]. Indeed, data are noisy and we choose not to give any prior knowledge on the input error structure, so this technique seems convenient. At this point, model of measurement and inversion have been developed. A question that naturally rises from here is how perturbations in the input space parameters affect the reconstruction of the size distribution. Figure 3 illustrates the necessity of sensitivity analysis for the all set of parameters.

Fig. 3. Sensitivity of the estimate of the aerosol size distribution when flow rate Q sh is varied from ±10%. The mode value is estimated with an absolute error around 6%.

IV. Uncertainty analysis In this part, a first uncertainty analysis is carried out to determine a confidence region on a simulated aerosol size distribution. Let's consider M : R d → R J the modelling function such that:

N = M(θ, R),
where N ∈ R J is the estimated aerosol size distribution vector and θ = (θ 1 , θ 2 , ..., θ d ) is the d-dimensional real vector of input parameters (geometrical parameters, environmental parameters, physical parameters, etc.) that is then decomposed in two sets of parameters: uncertain parameters denoted as vector θ u and fixed parameters with vector θ f . So, θ becomes θ = (θ f , θ u ) with

θ f = (θ f,1 , ..., θ f,d f ) ∈ R d f and θ u = (θ u,1 , ..., θ u,du ) ∈ R du , d = d u +d f .
As a first approach, uniform probability distributions are assigned to uncertain parameters θ u to model their associated uncertainty as follows:

θ u,i ∼ U([a i , b i ]) , i = 1, ..., d u .
However, the number of input parameters is high (d ∼ 40), so comes the priority to reduce the problem's dimension. It can be done with global sensitivity analysis methods. Since n(d p ) is not a scalar but a function of the mobility diameter, we pick some convenient functional as peaks width or peaks position. They play the role of surrogate measurements used for the global sensitivity analysis. According to the screening method called OAT Morris method [START_REF] Morris | Factorial sampling plans for preliminary computational experiments[END_REF], parameters that have been found to impact the most the reconstruction of the size distribution are denoted as the vector θ m u , θ m u = (Q sh , Q a , Q cpc , Q exc , r 1 , r 2 , L DM A , ∆ t , γ). Since the dimension of the problem has been significantly reduced, Monte Carlo simulations can now be performed. Now that significant parameters have been discriminated, others are fixed. Let θ f u these part of uncertain parameters considered as fixed, then new vector θ u can be written as θ u = (θ f u , θ m u ). If exponent [k] denotes the k-th Monte Carlo simulation, then

N [k] = M(θ f , (θ f u , θ m,[k] u ), R).
Figure 4 represents a 95% confidence interval based on uncertainty propagation results when an atmospheric aerosol is simulated. According to our knowledge, this is the first time that a confidence region is associated to an aerosol size distribution. Monte Carlo simulations have become possible in this case thanks to requirements expressed in terms of modeling and inversion.

V. Discussion A SMPS measurement model and associated uncertainties has been presented. They allow SMPS users to check their inversion method and also provide them with an idea of how uncertain the estimation of the aerosol size distribution can be. The uncertainty propagation is performed through the full model (measurement and inversion). However, In order to answer more precisely to SMPS users' need, the uncertainty analysis could be rethought as a decomposition of the global model into two sub-models. Moreover, improvements have to be brought concerning the probability distributions assigned to every uncertain parameter in order to get closer to the physical reality. Following this statement, a counting noise will also be added to the raw data response function.

Fig. 1 .

 1 Fig. 1. Scheme of the simulation principle

Fig. 2 .

 2 Fig. 2. Schematic diagram of a cylindrical DMA column

Fig. 4 .

 4 Fig. 4. 95% Confidence region on an atmospheric aerosol size distribution

  charge correction is used afterwards to estimate aerosol size distribution. Indeed, since electrical mobility function described in (2) is not bijective:∃(p , d p ) such that Z d (p , d p ) = Z d (p, d p ), large particles carrying multiple charges have the same mobility as smaller particles carrying less charges. n(d p ) is to describe aerosol particles carrying only one charge, so is the necessity of the correction.