Source Separation in Chemical Analysis:

Recent Achievements and Perspectives Leonardo T. Duarte Member, IEEE,, Saïd Moussaoui, Member, IEEE, and Christian Jutten, Fellow, IEEE

I. INTRODUCTION

Since its origins in the mid-1980s, the field of blind source separation (BSS) has attracted considerable attention within the signal processing community. One of the main reasons for such popularity is the existence of many problems that can be addressed in a BSS framework. Two noteworthy examples of applications can be found in audio and biomedical signal processing, for which a number of efficient solutions are now available. Additionally, there are relevant BSS problems in other domains but for which less effort has been put on. In this article, we deal with one of these fields, namely the field of analytical chemistry (AC), whose goal of is to identify or quantify, or both, chemical components present in a given analyte, i.e. the sample under analysis. As recently discussed in [START_REF] Duval | An overview of signal processing issues in chemical sensing[END_REF], several tasks in AC keep some relationship to the broad classes of detection and estimation problems typically found in signal processing.

Source separation is one of the most relevant estimation problems found in chemistry. Indeed, dealing with mixtures is paramount in different kinds of chemical analysis. For instance, there are some cases where the analyte is a chemical mixture of different components, e.g., in the analysis of rocks and heterogeneous materials through spectroscopy. Moreover, a mixing process can also take place even when the components are not chemically mixed. For instance, in ionic analysis of liquid samples, the ions are not chemically connected, but, due to the lack of selectivity of the chemical sensors, the acquired responses may be influenced by ions that are not the desired ones. Finally, there are some situations where the pure components cannot be isolated chemically since they appear only in the presence of other components. In this case, BSS may provide these components that cannot be retrieved otherwise.

In AC, there is a trend that advocates acquiring data through sensor arrays or by other instruments that are able to exploit diversity, which is essential in the application of source separation methods. Especially, Chemometrics [START_REF] Miller | Statistics and Chemometrics for Analytical Chemistry[END_REF], a sub field of AC, aims at extracting the relevant information from the multivariate signals provided by the chemical measurement sensor. Chemometrics is an established domain of research based on different backgrounds (physics, statistics, computer science, etc.). Several approaches that are now popular in signal processing arose in chemometrics -a typical example is a problem introduced in the 1970's known as multivariate curve resolution (MCR) [START_REF] Lawton | Self-modeling curve resolution[END_REF], which bears strong resemblance to non-negative matrix factorization (NMF) [START_REF] Lee | Learning the parts of objects by non-negative matrix factorization[END_REF].

Classically, the available methods in chemometrics for dealing with mixtures work in a supervised fashion, thus requiring a set of training samples. The application of supervised methods in this context has been proving successful in tasks such as odor and taste automatic recognition systems (electronic noses and tongues, respectively). However, despite this success, this approach suffers from, at least, two important practical problems. Firstly, the acquisition of training samples is usually a cost and timedemanding task. Secondly, due to the drift in the response of chemical sensors, which is usually caused by aging and variation of the acquisition conditions (temperature, pH, etc.), the calibration procedure must be performed before each new chemical analysis. Due to these limitations, there is a growing interest in solutions in which the calibration stage may be eliminated or, at least, simplified. This can be achieved by methods, like source separation or multi-way data analysis, able to exploit the measurement diversity for discovering latent variables. The goal in considering these approaches is to ultimately have chemical sensing systems that work in a plug-and-measure fashion, which may pave the way for devices operating in more challenging situations -for example, a calibration-free system would be quite helpful in the context of invasive physiological monitoring by means of miniaturized sensors.

Given the panorama described above, our aim is to shed some light on the use of BSS in chemical analysis. In this context, we firstly provide a brief overview on source separation (Section II), with particular attention to the classes of linear and nonlinear mixing models (Sections III and IV, respectively).

Then, (in Section V), we will give some conclusions and focus on challenging aspects that are found in chemical analysis. Although dealing with a relatively new field of applications, this article is not an exhaustive survey of source separation methods and algorithms, since there are solutions originated in closely related domains (e.g. remote sensing and hyperspectral imaging) that suit well several problems found in chemical analysis. Moreover, we do not discuss the supervised source separation methods, which are basically multivariate regression techniques, that one can find in chemometrics.

(n) = [x 1 (n), x 2 (n), . . . , x M (n)]
T , are given by:

x(n) = A(s(n)), (1)
where A(•) represents the operator associated to the mixing process, and the vector s 1) Source assumptions: Source separation can be solved with independent component analysis (ICA) [START_REF] Jutten | Blind separation of sources, part i: An adaptive algorithm based on neuromimetic architecture[END_REF][START_REF] Comon | Independent component analysis, a new concept?[END_REF], in which the sources are modeled as mutually independent random variables. Under some mild conditions, including the non-gaussianity of the sources, the application of ICA ensures a correct retrieval of the sources. It is also possible to achieve source separation by exploiting properties other than statistical independence, e.g. sparsity or non-negativity [START_REF] Comon | Handbook of blind source separation: independent component analysis and applications[END_REF]. Moreover, one can separate sources that are only mutually uncorrelated provided that additional temporal priors on the source signals hold, e.g., nonstationary or temporally correlated sources. This is very important in chemical applications since there are some cases in which the source independence assumption fails. For instance, when a chemical reaction takes place, the components exhibit dependency since either they are made of the same chemical elements or their concentrations vary according to stoichiometric coefficients, which stands for the coefficients in a balanced chemical equation.

(n) = [s 1 (n), . . . , s N (n)]
In chemical analysis, the sources present some features that can be exploited when developing a novel separation method. For instance, in ionic activity analysis (the effective concentration of an ion), the source signals correspond to time-series associated with the activities of each ion within the solution under analysis. Therefore, a first remarkable aspect is that the sources here are always non-negative, since there is no physical meaning in having a negative concentration. A second prior information is the fact that source signals are smooth, i.e. successive samples are temporally correlated. Indeed, concentrations usually present slow variations of amplitude along time. In order to illustrate these features, we show in Another category of chemical source signals are those obtained in optical spectroscopy (Raman, Infrared, etc.). Again, there are some valuable information that can be exploited: (1) the sources are again non-negative, since they correspond to absorption, reflectance or diffusion spectra, (2) there exist many libraries containing the spectra of different chemical components. This information can be used to adapt BSS algorithms based on sparse representation over a huge dictionary whose atoms are the pure component spectra. However, it would be difficult to use a reference library since the component spectra may present some variability in actual experiments depending on parameters such as temperature or pH condition. In Figure 1(b), an example of spectral source (infrared reflectance of calcite) is provided.

2) Mixing models: A lot of effort was undertaken to develop methods tailored to linear and instantaneous mixing systems in which the number of sources is equal to the number of mixtures (N = M ).

Practically, this assumption implies that the number of sources is known, which is often a tricky issue.

In this case, the operator A(•) is given by a square matrix A ∈ R N ×N . Moreover, there are many works that deal with the convolutive case, in which each entry of the mixing matrix corresponds to a linear filter. Finally, there is an important topic in BSS where the aim is to set up separation methods in the case of nonlinear mixing models. In this context, one must bypass many difficulties that do not arise in the linear case. For instance, ICA does not ensure separation in a general nonlinear context, that is, retrieving independent components is not enough to reconstruct the actual sources. As a consequence, December 19, 2013 DRAFT some researchers have been considering constrained classes of nonlinear model for which ICA is still valid. For instance, ICA-based solutions ensure source separation in an important class of constrained models known as post-nonlinear (PNL) models ( [START_REF] Comon | Handbook of blind source separation: independent component analysis and applications[END_REF], chapter 14). In chemical mixture analysis, one can find in the literature works on both linear and nonlinear mixing models. In the sequel, we provide an (non-exhaustive) overview of source separation problems arising in each situation.

III. LINEAR MODELS

Measurement techniques such as optical spectroscopy and/or chromatography are frequently used in chemical analysis to extract relevant information related to the chemical composition of an heterogeneous material (identify the components and assess their relative abundance [START_REF] Malinowski | Factor Analysis in Chemistry[END_REF]). The data processing corresponds to a source separation problem where the linear instantaneous mixture model holds with relatively small concentrations of the components in absorption spectroscopy, thanks to the Beer-Lamber law [START_REF] Pfeiffer | The origins of Beer's law[END_REF].

In such model, each measured spectrum x i (λ) is a linear combination of the component spectra s j (λ),

and the mixing coefficients a ij are related to the abundances of each component:

x i (λ) = N ∑ j=1 a ij s j (λ). (2)
In practice, M spectra (i.e., i = 1 . . . M in (2)) are measured for different chemical conditions and L values of a spectral variable (frequency, wavelength or wavenumber) λ. In order to develop a linear spectral mixture separation method, an objective function must be obtained by formulating some assumptions on the sought component spectra (the sources) and on their abundances (the mixing coefficients). Then, one must define a mathematical algorithm to optimize this objective function. The main information in spectral mixture data analysis is the non-negativity of both the pure component spectra and abundances. However, even with a linear mixing model, only accounting for these constraints does not lead to a unique solution [START_REF] Lawton | Self-modeling curve resolution[END_REF]. Additional constraints or assumptions are thus used to select a particular (and meaningful) solution among the feasible ones [START_REF] Tauler | Simultaneous analysis of several spectroscopic titrations with selfmodelling curve resolution[END_REF]. The existing methods differ on the way the constraints and the additional assumptions are formulated and how the separation is performed. A short overview of the different approaches is firstly given in this section and an example illustrates the application of one of these techniques to the separation problem.

A. Independent Component Analysis

As mentioned before, since BSS is an ill-posed inverse problem, prior knowledge and/or additional assumptions should be used to get a unique and correct solution. Principal component analysis (PCA) [START_REF] Hotelling | Analysis of a complex of statistical variables into principal components[END_REF] assumes that the signals to be reconstructed are mutually uncorrelated but this orthogonality constraint does not ensure the non-negativity of the solution. A more constraining statistical assumption used for source separation is the mutual independence of the sources, leading to the ICA concept [START_REF] Comon | Handbook of blind source separation: independent component analysis and applications[END_REF], for which many algorithms have been developed (see [START_REF] Comon | Handbook of blind source separation: independent component analysis and applications[END_REF]). If the sought source signals are mutually statistically independent, they can be separated successfully by ICA methods and their non-negativity will be ensured implicitly (at least, only few negative values appear in the estimates) as reported in [START_REF] Nuzillard | Application of blind source separation to 1-D and 2-D Nuclear Magnetic Resonance spectroscopy[END_REF], where a December 19, 2013 DRAFT second-order blind identification (SOBI) algorithm [START_REF] Belouchrani | A blind source separation technique using second-order statistics[END_REF] was applied to the analysis of nuclear magnetic resonance (NMR) data. But, when the source signals are not mutually independent or when their mutual independence is not observed due to the finite and small number of samples, the non-negativity information should be considered. In [START_REF] Plumbley | Algorithms for non-negative independent component analysis[END_REF], Plumbley showed that it is possible to incorporate jointly non-negativity and mutual independence of the sources (See also in [START_REF] Comon | Handbook of blind source separation: independent component analysis and applications[END_REF], chapter 13). This method yields a correct solution providing the condition of well-grounded sources (that is sources having non-vanishing probability distributions around zero) is fulfilled.

B. Multivariate Curve Resolution

This approach, proposed by Lawton and Sylvestre [START_REF] Lawton | Self-modeling curve resolution[END_REF], termed by multivariate curve resolution, firstly decomposes the data matrix, using Principal Component Analysis (PCA), X = U V . Then, a linear transformation T is calculated in order to transform the principal components V and their weight matrix U into non-negative estimates of pure spectra S = T V and mixing coefficients A = U T -1 , respectively.

Since accounting for non-negativity alone does not ensure the uniqueness of the solution, this approach leads to several possible values of matrix T which provide the set of admissible (feasible) solutions [START_REF] Lawton | Self-modeling curve resolution[END_REF].

In order to reduce this set Sasaki et al. [START_REF] Sasaki | Constrained nonlinear method for estimating component spectra from multicomponent mixtures[END_REF] suggests to add further constraints, in addition to the nonnegativity, and proposes to search a linear transformation T by minimizing a two-term criterion, in which a first part penalizes negative estimates of the pure spectra and mixing coefficients, and the second part uses an entropic cost function to make the estimated spectra smoother and mutually independent. The resulting optimization is given by:

min T   ∥ min(S, 0)∥ 2 F + ∥ min(A, 0)∥ 2 F + β L ∑ n=1 N ∑ j=1 s j (n) log s j (n)   ,
where ∥ • ∥ F stands for the Frobenius norm and β is a regularization parameter that allows to adjust the trade-off between the two parts. The optimization of the whole objective function is performed by using the Nelder-Mead algorithm. However, this method may converge to local or spurious minima since this criterion in non-convex and the criterion shape highly depends on the regularization parameters that should be specified manually. This method was revisited in [START_REF] Widjaja | Pure component spectral reconstruction from mixture data using SVD, global entropy minimization, and simulated annealing. numerical investigations of admissible objective functions using a synthetic 7species data set[END_REF] where a simulated annealing optimization algorithm is used. In the signal processing community, the algebraic approach for ICA is very similar to that of Lawton-Sylvestre but it is based on the statistical independence of the source signals (pure spectra). To get mutually independent signals, the transformation matrix is reduced to a unitary rotation matrix by minimizing a contrast function ( [START_REF] Comon | Handbook of blind source separation: independent component analysis and applications[END_REF], chapter 3). Assuming the independence of the sources leads to a unique solution but it does not guarantee its non-negativity.

December 19, 2013 DRAFT

C. Non-negative Matrix Factorization

This approach performs a constrained least squares estimation of S and A

( Ŝ, Â) = argmin S,A ||X -AS|| 2 F .
This (non-convex) optimization problem can be solved using alternating least square (ALS) with alternate exact solving [START_REF] Tauler | Simultaneous analysis of several spectroscopic titrations with selfmodelling curve resolution[END_REF] or multiplicative [START_REF] Lee | Learning the parts of objects by non-negative matrix factorization[END_REF] alternate update methods. At each iteration k, ALS algorithms minimize alternatively the above criterion with respect to S, keeping A fixed, or to A, keeping S fixed, with non-negativity constraints [START_REF] Bro | PARAFAC. tutorial and applications[END_REF] on S and A. It leads to solve the following optimization problems

Ŝ(k) = argmin S ||X -Â(k-1) S|| 2 F s.t. S 0 (3)
Â(k) = argmin A ||X -ÂS (k) || 2 F s.t. A 0. (4)
Multiplicative methods update iteratively the estimates of the sources and the mixing coefficients using a particular multiplicative learning rule that ensures the non-negativity [START_REF] Lee | Learning the parts of objects by non-negative matrix factorization[END_REF]. However since the criterion is non-convex, NMF algorithms do not lead a unqiue solution, unless in some particular conditions [START_REF] Donoho | When does non-negative matrix factorization give a correct decomposition into parts?[END_REF].

Therefore, the NMF results highly depend on the algorithm initialization. Several contributions proposed to initialize the algorithm by the results obtained with an unconstrained decomposition method such as PCA, factor analysis algorithms, or using pure variable detection methods such as simple-to-use interactive self modeling mixture analysis (SIMPLISMA) [START_REF] Windig | Interactive self-modeling mixture analysis[END_REF] and orthogonal projection approach (OPA) [START_REF] Sanchez | Orthogonal projection approach applied to peak purity assessment[END_REF]. In [START_REF] Nuzillard | Application of blind source separation to 1-D and 2-D Nuclear Magnetic Resonance spectroscopy[END_REF], for NMR spectroscopy, it has been shown that ICA methods can be used successfully for the initialization of the ALS approach. However, each initialization leads to a local minimum of the criterion and none of these methods has proven to outperform the others since the result of each method highly depends on the data at hand.

Similarly to the case of curve resolution methods, additional constraints such as sum-to-one (also called closure or full additivity) and unimodality (i.e. presence of only one maximum in each column of matrix A) may be added to reduce the set of admissible solutions [START_REF] Tauler | Calculation of maximum and minimum band boundaries of feasible solutions for species profiles obtained by multivariate curve resolution[END_REF]. The sum-to-one constraint corresponds to assuming that the sum of the elements of each row of matrix A is equal to one (or to a known constant). This is the case for instance in reaction-based systems, where a mass balance equation is obeyed by the concentration profiles of the species present in the system. See for instance [START_REF] De Juan | Assessment of new constraints applied to the alternating least squares method[END_REF] 1 . Another constraint that reduces the set of admissible solution is the local rank or selectivity [START_REF] Tauler | Calculation of maximum and minimum band boundaries of feasible solutions for species profiles obtained by multivariate curve resolution[END_REF] which refers to the fact that for certain rows or columns of the data matrix, some mixing coefficients are known to be non-zero while other coefficients are known to have zero values.

There are other constraints that can be taken into account through a penalized least squares estimation approach. This idea is used in positive matrix factorization (PMF) [START_REF] Paatero | Positive matrix factorization: a nonnegative factor model with optimal utilization of error estimates of data values[END_REF] (this work was historically the first one dealing with spectral mixture analysis as a matrix factorization problem under positivity constraints)

where the criterion to minimize is min

S,A   ||X -AS|| 2 F + α||A|| 2 F + β||S|| 2 F -γ L ∑ n=1 N ∑ j=1 log s j (n) -δ M ∑ i=1 N ∑ j=1 log a ij   . (5)
The hyperparameters γ and δ control the strength of the logarithmic barrier function that prevents negative values of the source samples and mixing coefficients. Hyperparameters α and β allows to adjust the weight of the quadratic regularization criteria. Contrarily to the constrained least squares methods, the PMF approach leads to an unconstrained optimization problem.

It is also possible to consider the sparsity constraint through a penalized least squares approach. Such an idea gives rise to non-negative sparse coding (NNSC) [START_REF] Hoyer | Non-negative matrix factorization with sparseness constraints[END_REF], which searches for minimizing the following

This approach was applied by [START_REF] Sajda | Nonnegative matrix factorization for rapid recovery of constituent spectra in magnetic resonance chemical shift imaging of the brain[END_REF] for the analysis of magnetic resonance chemical shift imaging data.

There are other alternative approaches for reducing the set of admissible solutions, see for instance [START_REF] De Juan | Local rank exploratory analysis of evolving rank-deficient systems[END_REF] and [START_REF] Miao | Endmember extraction from highly mixed data using minimum volume constrained nonnegative matrix factorization[END_REF].

D. Bayesian Approaches

The formulation of source separation using Bayesian estimation theory is reported in details in ( [START_REF] Comon | Handbook of blind source separation: independent component analysis and applications[END_REF],

chapter 12). Its first application to the separation of spectral mixture data has been proposed in [START_REF] Moussaoui | Separation of non-negative mixture of non-negative sources using a Bayesian approach and MCMC sampling[END_REF]. The where N (µ, Σ) stands for the Gaussian distribution with mean µ and covariance matrix Σ. This matrix reduces to σ 2 I L when the noise samples are assumed to be independent and identically distributed. In order to account for the non-negativity and the sparsity of the component spectra, a Gamma distribution model was used in [START_REF] Moussaoui | Separation of non-negative mixture of non-negative sources using a Bayesian approach and MCMC sampling[END_REF], while a truncated Gaussian distribution or a Dirichlet distribution can be used to encode the sum-to-one constraint on the mixing coefficients [START_REF] Dobigeon | Semi-supervised linear spectral unmixing using a hierarchical Bayesian model for hyperspectral imagery[END_REF].

The key point of the Bayesian approach is to apply Bayes' theorem to express the a posteriori density,

p(S, A|X, θ, Σ) = p(X|A, S, θ, Σ) p(A|θ) p(S|θ, Σ)/p(X), (7)
from which a statistical inference is conducted to perform the separation. This distribution combines explicitly the available assumptions on the pure component spectra and their abundances with the information coming from the measured data. The inference of the unknown quantities can be conducted by minimizing J(S, A) =log p(S, A|X, θ) or calculating empirical posterior means of S and A after drawing samples from p(S, A|X, θ) using Monte Carlo Markov Chain (MCMC) methods (See [START_REF] Moussaoui | Separation of non-negative mixture of non-negative sources using a Bayesian approach and MCMC sampling[END_REF][START_REF] Dobigeon | Semi-supervised linear spectral unmixing using a hierarchical Bayesian model for hyperspectral imagery[END_REF] for details). The first approach is equivalent to minimizing a criterion similar to (6) with regularization criteria R(A) and R(S) linked to the prior distributions according to R(A) =log p(A|θ) and R(S) =log p(S|θ), while the second approach allows to infer all the weighting parameters at the price of a significant increase of the computation time. The Bayesian method was successfully applied for the analysis of Hyperspectral data in [START_REF] Moussaoui | On the decomposition of Mars hyperspectral data by ICA and Bayesian positive source separation[END_REF] and chemical reaction monitoring [START_REF] Carteret | Polymorphism studied by lattice phonon Raman spectroscopy and statistical mixture analysis method[END_REF].

E. Geometrical Methods

Geometrical methods are based on the empirical distribution (geometrically speaking, the scatter plot) of the mixture data. Since these data result from a non-negative mixing of non-negative data, the scatter plot of the mixed data is contained in the simplicial cone generated by the columns of the mixing matrix [START_REF] Donoho | When does non-negative matrix factorization give a correct decomposition into parts?[END_REF]. Figure 3 illustrates three examples of mixture data distribution for different types of sources. It can be noted that in the case of sparse sources, each row of S has a dominant peak at some location (column number) where other rows have zero elements, then the problem of finding the columns of the mixing matrix A reduces to the identification of the edges of a similicial cone, edges where the data are concentrated 3(b). In the case of sum-to-one mixing, the mixture data will be distributed on a simplex whose vertices correspond to mixture data containing only one component. component. Thus, each vertice is associated with a column of the mixing matrix A. In [START_REF] Ouedraogo | Simplicial cone shrinking algorithm for unmixing nonnegative sources[END_REF], efficient algorithms are designed for data where the sources are not well grounded and the mixing matrix does not satisfy the sum-to-one or the December 19, 2013 DRAFT pure pixel assumptions. This method can also handle noisy data and was applied to Positron emission tomography (PET) imaging and mass spectroscopy. In chemical imaging spectrometry and remote sensing by hyperspectral imaging [START_REF] Chang | Hyperspectral Data Exploitation[END_REF][START_REF] Bioucas-Dias | Hyperspectral unmixing overview: Geometrical, statistical, and sparse regression-based approaches[END_REF], the spectral mixture analysis is often handled using a two step procedure: the pure spectra estimation and the abundance fraction assessment, respectively. In the first step, the pure component of the mixture are identified by using an Endmember Extraction Algorithm (EEA). See for instance [START_REF] Plaza | A quantitative and comparative analysis of endmember extraction algorithms from hyperspectral data[END_REF] for a recent performance comparison and discussion on EEAs. The most popular EEAs is the N-FINDR algorithm.

N-FINDR estimates the pure component spectra by identifying the largest simplex whose vertices are taken from the convex hull of the data. Another popular and faster alternative is the (Vertex Component Analysis-VCA) method which has been proposed in [START_REF] Nascimento | Vertex component analysis: a fast algorithm to unmix hyperspectral data[END_REF]. It consists in iteratively estimating the vertices of the simplex without calculating the convex hull. A common assumption in VCA and N-FINDR is the existence of pure pixels (pixels composed of a single component) in the observed data. Alternatively, the Minimum Volume Transform (MVT) finds the smallest simplex that contains all the pixels [START_REF] Craig | Minimum-volume transforms for remotely sensed data[END_REF][START_REF] Miao | Endmember extraction from highly mixed data using minimum volume constrained nonnegative matrix factorization[END_REF]. The second step of the spectral unmixing can use various strategies such as those based on constrained least squares estimation [START_REF] Moussaoui | Primal-dual interior point optimization for penalized least squares estimation of abundance maps in hyperspectral imaging[END_REF]. one or more additional diversities can be considered, leading to a 3-way (or more generally a multi-way) array of data, which can be represented by a tensor. As an example, in fluorescence spectroscopy, the (measured) fluorescence intensity dependents on 3 variables: the fluorescence emission spectrum, the absorbance spectrum and the relative concentration of the components. Of course, tensorial decomposition is not restricted to fluorescence, and other applications include time resolved spectroscopy [START_REF] Niazi | Parafac decomposition of three-way kinetic-spectrophotometric spectral matrices based on phosphomolymbdenum blue complex chemistry for nitrite determination in water and meat samples[END_REF],

F. Tensorial methods

Most

multidimentional NMR [START_REF] Nilsson | T1-diffusion-ordered spectroscopy: Nuclear magnetic resonance mixture analysis using parallel factor analysis[END_REF] or polarized Raman spectroscopy [START_REF] Miron | Joint processing of the parallel and crossed polarized raman spectra and uniqueness in blind nonnegative source separation[END_REF].

Historically, PARAFAC or canonical polyadic decomposition (CPD) -inspired from early works on factor analysis in psychometrics -has been intensively studied in chemometrics from 90's and popularized by R. Bro [START_REF] Bro | PARAFAC. tutorial and applications[END_REF]. Currently, many theoretical contributions and applications are addressing CPD decomposition of 3-way arrays. More details and references of tensorial methods and an application in fluorescence spectroscopy are provided, respectively, in the articles [START_REF] Zhou | Algorithms for nonnegative matrix and tensor factorizations of large-scale and multi-block data[END_REF] and [START_REF] Comon | Tensor decompositions: a brief survey[END_REF], found in the present special issue.

G. Application Example

Calcium carbonate is a chemical material used for a large variety of applications such as filler for Moreover for each temperature, Raman spectra were collected at regular time intervals for monitoring the phase transformation. Finally, a total of M = 37 spectra (each one obtained for a given temperature and phase transformation time) of L = 477 wavelengths are obtained. Details on the experiment can be found in [START_REF] Carteret | Polymorphism studied by lattice phonon Raman spectroscopy and statistical mixture analysis method[END_REF].

IV. NONLINEAR MODELS

Although linear mixing models provide a satisfactory first-order approximation of the measurement process in several applications, there are some cases in which a nonlinear model is mandatory. Indeed, exploiting the physical theory behind the measurement transducer allows to increase the amount of information that can be extracted from the data. For example, the sensitive membrane of potentiometric sensors such as ion-selective electrodes (ISEs) or ion-selective field effect transistors (ISFETs) can be described by the classical Nernst equation [START_REF] Gründler | Chemical sensors: an introduction for scientists and engineers[END_REF], which establishes a logarithmic relation between ionic activity and the membrane potential. In this section we shall discuss how source separation methods can be set to deal with nonlinear mixing models, paying special attention to the illustrative example of ISE arrays. Afterward, some considerations on other types of sensors are made.

ISEs, which are the most used chemical sensors in industry, are simple devices for measuring ionic activity. A well know example of ISE is the glass electrode, which is used for measuring the pH of a given solution. Besides, one can find ISEs tailored to different ions such as ammonium, potassium and sodium. These devices have been intensively used, for instance, in food and soil inspection, clinical analysis and water quality monitoring. One of the reasons that explains the success of ISEs in such applications is the simplicity of this approach. Indeed, ISE-based analyses do not require sophisticated laboratory equipments and procedures and, thus, can be carried out in the field if necessary. Although attractive, potential electrodes such as ISEs and ISFETs suffer from an important drawback: they are not selective, that is, the generated potential depends on a given target ion but also on other undesirable ions, which are called interfering ions [START_REF] Gründler | Chemical sensors: an introduction for scientists and engineers[END_REF]. There are some situations in which this interfering process is so weak that it can be neglected. However, this phenomenon may become important when the target ion and the interfering ones have similar physical and/or chemical properties. A possible solution to eliminate this interference effect is to consider the diversity provided by an ISE array using a source separation approach.

The interference in ISE can be modeled through an empirical extension of the Nernst equation called the Nicolsky-Eisenman (NE) equation [START_REF] Gründler | Chemical sensors: an introduction for scientists and engineers[END_REF]. According to this model, the response of the i-th sensor December 19, 2013 DRAFT (dedicated to measure the i-th ion) within the array is given by

x i (n) = e i + d i log   s i (n) + ∑ j,j̸ =i a ij s j (n) zi/zj   , (8)
where d i and e i are constants that depend on some physical parameters, z i and s i (n) denote, the valence and the activity of the i-th ion, respectively. The non-negative parameters a ij , the selectivity coefficients, explain the influence (interference) of the j-th ion on the i-th sensor; n stands for the temporal index.

There are some pairs of ions for which high values of selectivity coefficients are observed. When the valences of the ions are different, the model ( 8) becomes difficult to deal with because a nonlinearity (power term) appears inside the logarithm term [START_REF] Duarte | A nonlinear source separation approach to the Nicolsky-Eisenman model[END_REF]. However, when the valences are equal, which is the most common situation, Equation ( 8) becomes a special nonlinear (NL) model known as post-nonlinear (PNL) model [START_REF] Duarte | A Bayesian nonlinear source separation method for smart ion-selective electrode arrays[END_REF]. Indeed, in PNL systems, there is a linear mixing stage followed by a set of componentwise nonlinear function, which in the case of ISE arrays correspond to logarithms. As mentioned before, an interesting property of PNL models is that, under conditions that are close to those required in the linear case [START_REF] Comon | Handbook of blind source separation: independent component analysis and applications[END_REF], the application of ICA leads to source separation.

A. Nonlinear Independent Component Analysis

For equal valences, when the number of mixtures and sources is the same, the mixing process described in [START_REF] Duarte | A dataset for the design of smart ion-selective electrode arrays for quantitative analysis[END_REF] can be counterbalanced by a separating system in which each estimated source is given by

y k (n) = M ∑ i=1 w ki exp(f i + h i x i (n)), for k = 1, . . . , N (9)
where w ki , f i and h i are unknown parameters that must be adjusted and M the number of sensors and N the number of target ions-actually, the parameters f i only introduce a scale ambiguity, which is usual in BSS, and, therefore, can be fixed a priori to a given value.

As in the linear case, the adjustment of PNL separating models can be performed by means of ICA, but some care must be taken while developing the separation algorithm. For instance, one must resort to stronger independence measures such as the mutual information. The adjustment of the separating system is thus carried out by the following optimization problem min wki,hi

I(y 1 (n), . . . , y N (n)), (10)
where I(y 1 , . . . , y N ) denotes the mutual information between the estimated signals y(n). There are several algorithms for solving [START_REF] Pfeiffer | The origins of Beer's law[END_REF], where main difficulties are related to the estimation of the mutual information and to the existence of spurious local minima in this cost function. In [START_REF] Duarte | Blind source separation of postnonlinear mixtures using evolutionary computation and order statistics[END_REF], for instance, December 19, 2013 DRAFT the optimization problem expressed in [START_REF] Pfeiffer | The origins of Beer's law[END_REF] was tackled by a bio-inspired optimization method that is robust to convergence to local minima and that does not require the calculation of the derivatives of [START_REF] Pfeiffer | The origins of Beer's law[END_REF].

Moreover, [START_REF] Duarte | Blind source separation of postnonlinear mixtures using evolutionary computation and order statistics[END_REF] adopts a mutual information estimator based on order statistics.

In order to illustrate the application to a real situation of the ICA method [START_REF] Duarte | Blind source separation of postnonlinear mixtures using evolutionary computation and order statistics[END_REF], we considered an experiment where the data was acquired by an array composed of a sodium and a potassium ISE (two mixtures). Therefore, there are two mixtures and two sources (the activities of potassium and sodium).

Data acquisition was achieved via a flow-injection analysis system, which usually increases the stability and the sensitivity of the array. In Figure 5

B. Bayesian Approach

Bayesian source separation can also be applied in the context of ISE arrays. Actually, differently from ICA, whose basis lies in the independence assumption, the Bayesian approach does not introduce a measure of independence but it searches a solution allowing to jointly explain the data according to the mixing model and fulfill some prior knowledge. Similarly to the linear case (see Equation ( 7)), the posterior distribution in the nonlinear case is given by p(s, A, g|x) = p(x|A, s, g)p(A)p(s)p(g)/p(x),

where the vector g includes the parameters associated with the logarithmic functions (see Equation ( 8)).

As in the linear case, the inference step can be done based on a set of samples drawn from the posterior distribution via MCMC methods. The main difficulty when dealing with nonlinear models using MCMC methods is the need to draw samples from distributions for which classical sampling method may be inefficient. In [START_REF] Duarte | A Bayesian nonlinear source separation method for smart ion-selective electrode arrays[END_REF], by relying on lognormal prior distributions on p(s), limited support gaussian densities on p(g), we showed the relevance of accounting for the the non-negativity constraint and introducing auxiliary variables to handle the sampling difficulties.

C. Exploiting other prior information

In addition to the non-negativity of the sources, other prior information can be considered when developing nonlinear source separation methods. For example, in [START_REF] Duarte | Blind compensation of nonlinear distortions: Application to source separation of post-nonlinear mixtures[END_REF], a method exploiting the assumption that December 19, 2013 DRAFT the sources are bandlimited is proposed to compensate the nonlinear stage of PNL systems. The advantage brought by this approach is that, once the nonlinear component-wise functions are counterbalanced, the resulting separation problem becomes linear and, therefore, can be dealt with by means of linear BSS methods. In the context of ISE arrays, assuming bandlimited sources is realistic since the signals of interest have a spectral content concentrated on low-frequency bands. The approach proposed in [START_REF] Duarte | Blind compensation of nonlinear distortions: Application to source separation of post-nonlinear mixtures[END_REF] was applied in the analysis of the ions ammonium and potassium, leading to a good estimation of the sensors' nonlinearities.

Another prior information that can be considered is the existence of silent periods of a given ion, i.e., time windows in which the concentration of a given ion is approximately zero. Interestingly, this scenario is close to that of found in speech separation and can be interpreted as a sparsity-based approach. In [START_REF] Duarte | A nonlinear source separation approach to the Nicolsky-Eisenman model[END_REF], December 19, 2013 DRAFT a method that uses this prior was proposed to estimate the parameters d i of (8) in the case where the valences of the ions are different. In this situation, as can be noted in Equation ( 8), the mixing model is much more complex than PNL model, and, thus, estimating the nonlinearities at first usually simplifies the problem.

D. Linear-quadratic mixing model

In addition to the complex model ( 8) considered in the context of ISE arrays for analyzing ions with different valences, there is another relevant class of models that arises in chemical analysis: the linearquadratic (LQ) model ( [START_REF] Comon | Handbook of blind source separation: independent component analysis and applications[END_REF], chapter 14). An interesting aspect of the LQ mixture is that it can be seen as linear mixture of dependent sources. Such a feature simplifies the derivation of separation algorithms in this case. According to Clifford-Touma equation [START_REF] Yamazoe | Theory of power laws for semiconductor gas sensors[END_REF], the mixing process that takes place in the analysis of two gases by using a tin oxide gas sensor can be described by a LQ model. Moreover, LQ models also arise in the context of hyperspectral imaging, especially when there are multiple reflexions caused by the presence of buildings or trees [START_REF] Bioucas-Dias | Hyperspectral unmixing overview: Geometrical, statistical, and sparse regression-based approaches[END_REF]. This is an interesting problem since hyperspectral imaging is closely-related to the above-mentioned problems of separating spectra. Such a model is also considered in fluorescence spectroscopy of highly concentrated solution, involving a screen effect, and a possible method for solving this problem is proposed in [START_REF] Luciani | A simple correction method of inner filter effects affecting feem and its application to the parafac decomposition[END_REF].

V. CONCLUSIONS AND FURTHER CHALLENGES

In this article, we aimed at presenting an overview on how source separation methods can be applied in chemical analysis. We briefly discussed some kinds of chemical data, and presented several approaches that can be applied. A typical characteristic of chemical data is the existence of prior information on the sources and the mixing process. These priors have motivated researchers working on the subject to consider alternatives to the framework of ICA, thus showing that, besides being an exciting field of application to classical BSS methods, chemical analysis may also be inspiring environment for novel BSS paradigms. Another interesting aspect of this relatively new field of applications is related to the challenges that are found when developing separation methods, e.g. nonlinear mixtures and dependent sources.

Despite the current advances achieved in source separation for chemical data, there are still many issues that must be tackled.

1) Mixing models:

A first important point concerns the mixing models that should be adopted. Indeed, transducer physics gives realistic mixing models whose processing in a raw version or even after some December 19, 2013 DRAFT simplifications will provide very challenging signal processing problems. Moreover, most of the models give a static description, that is, we assume that the mixtures follow an instantaneous model. Therefore, if a dynamical mixing model were available, it would be possible to take advantage of the sensor's dynamics to obtain better separation results. This could be achieved by considering convolutive models.

2) Noise and measurement process modeling: Another difficult aspect that must be handled in chemical analysis is the existence of noise and other complex phenomena. In ISEs, for instance, since there is nonlinear elements in the transducer stage, there may be a strong noise amplification during the separation process -one can show that a separation method based on inversion leads to a multiplicative noise model [START_REF] Duarte | A Bayesian nonlinear source separation method for smart ion-selective electrode arrays[END_REF]. Moreover, in the context of spectral unmizing, finding realistic noise models for non-negative sources is still an open issue. Finally, when dealing with chemical sensors one may find hysteresis, which is a nonlinear and dynamical phenomenon that is very difficult to model.

3) Source number: the source number is one of the oldest questions in signal analysis and remains of huge interest for the application of source separation methods to real applications. In the linear mixing model, additive noise and mutually uncorrelated sources, it is generally addressed by a subspace analysis approach. However, the theory behind such analysis breaks down when dealing with either nonlinear mixing models or non-negative sources. For instance, by considering non-negativity in the linear mixing model the source number corresponds to the non-negative rank which is different from the classical matrix rank. Indeed, when dealing with mixture data related some chemical reaction monitoring, the number of sources will be known in advance.

4)

Enhancing the data acquisition setup: the performances of any signal processing method will highly depend on how much the data fulfill the method hypotheses. For instance, the application of NMF will be more efficient when the sources are well grounded, when there is some pure pixels in the spectral image and, more generally, when the uniqueness conditions are respected. It would be therefore useful to pay attention on how to design the measurement setup in such a way to meet, at best, theses conditions or how to introduce additional measurement modalities to make the mixture analysis problem more tractable.

In that respect, chemists play a fundamental role here especially in data acquisition. There is a great number of different acquisition techniques, even for the same kind of analysis, and many properties, such as temperature, pH, pressure, humidity, that when correctly controlled by a chemist, may lead to a good diversity within the acquired data. Therefore, in order to further advance in the application of source separation in chemistry, the collaboration between researchers in the fields of analytical chemistry and signal processing should be reinforced, even more than before.

December 19, 2013 DRAFT

Figure 1 (

 1 Figure 1(a) an example of sources related to the concentrations of potassium and ammonium ions and

  Infrared reflectance of calcite.

Fig. 1 .

 1 Fig. 1. Examples of sources in chemical analysis.

Figure 2 (

 2 b) illustrates three absorption spectra resulting from the linear mixing of two spectral sources, shown in Figure2(a). By considering a matrix notation for Equation (2) and a noise term, the measured spectra can be stored in the rows of the data matrix X ∈ R M ×L , that can thus be factorized according to X = AS + E. The pure component spectra are identified as the sources (rows of matrix S ∈ R N ×L ) and the abundance fraction of each component are the elements of the mixing matrix A ∈ R M ×N . Matrix E ∈ R M ×L represents measurement errors and any deviation from the linear mixing model. In the chemical analysis community, the problem of estimating S and A knowing only matrix X is termed by Multivariate Curve Resolution, Spectral Mixture Analysis and Factor Analysis while in signal and image processing field this problem is called Source Separation and Spectral Unmixing.

2 F

 2 + β||S|| 1 ) s.t. A 0 and S 0.

  spectral mixture separation in a Bayesian framework consist in describing the statistical properties of the measurement noise and assigning parametric a priori distributions p(S|θ) and p(A|θ) on the pure component spectra and abundances, respectively. These distribution are defined with some hyperparameters represented by θ. A multivariate Gaussian distribution is generally used to model the noise statistics, which yields the likelihood p(X|A, S, Σ) = N (X -AS; Σ),December 19, 2013 DRAFT

Fig. 3 .

 3 Fig. 3. Illustration of the simplicial cone and the empirical distribution of the mixture data for three types of non-negative sources : (a) uniform well grounded sources, (b) sparse sources and (c) uniform sources with sum-to-one mixing constraint.

  of the above methods simply exploit one type of diversity: it leads to a 2-way array of data, easily represented by a matrix. Above examples consider 2-way arrays based on space and time (mixtures of signals), or space and frequency (spectral mixtures) dimensions. But in many chemical experiments, December 19, 2013 DRAFT

  plastics or paper. Depending on operating conditions, Calcium carbonate crystallizes as Calcite, Aragonite or Vaterite. The formation of Calcium carbonate by mixing two solutions containing respectively Calcium and Carbonate ions takes place in two steps. The first step is the precipitation one, which is very fast and provides a mixture of Calcium carbonate polymorphs. The second step (a slow process) represents the phase transformation from the unstable polymorphs to the stable one (Calcite). The physical properties of the crystallised product depend largely on the polymorphic composition, so it is necessary to quantify these polymorphs when they are mixed. The main purpose of this application is to study the relation of polymorphs and temperature and to explore favorable conditions for Calcite formation.Calcium chloride and Sodium carbonate separately dissolved in Sodium chloride solutions of the same concentration were rapidly mixed to precipitate Calcium carbonate. A sample was collected 2 minutes after the beginning of the experiment to determine the polymorphic composition at the end of the precipitation step. Raman spectra of this sample have been measured for various temperatures ranging between 20 o C and 70 o C to determine the influence of temperature on the polymorphs precipitation.

Figure 4 (Fig. 4 .

 44 Fig. 4. Separation of Raman spectra of calcium carbonate: mixture data and separation results. (a) Collected mixture spectra at 2 minutes after the beginning of the phase transformation for different temperature values. (b) Comparison between the estimated sources (continuous) and the reference spectra (dashed) of the three components. (c) Abundances of the three components at the beginning of the phase transformation for different temperature values. (d) Temporal evolution of the three component abundances at T=40 o C.

Figure 4 (

 4 Figure 4(b) shows the estimated sources using a Bayesian separation approach with a Gamma distribution prior on the sources and a Dirichlet distribution as prior on the mixing coefficients[START_REF] Moussaoui | Separation of non-negative mixture of non-negative sources using a Bayesian approach and MCMC sampling[END_REF][START_REF] Dobigeon | Semi-supervised linear spectral unmixing using a hierarchical Bayesian model for hyperspectral imagery[END_REF]. From the spectroscopic point of view and according to the locations of the vibrational peaks, the identification of the three components is easy. The evolution of the concentration profile versus the temperature is shown in Figure4(c). It can be deduced that pure Vaterite is observed at 20 o C and a quite pure aragonite is obtained at 60 o C . However, between 20 o C and 60 o C ternary mixtures are observed. The abundance of calcite is maximal at 40 o C which is in agreement with results reported in the literature. Let us now

  (a), the actual sources are depicted. As a result of the lowselectivity of each sensor within the array, the acquired signals correspond to nonlinear mixtures of the sources, which is illustrated in Figure5(b). After solving the problem[START_REF] Pfeiffer | The origins of Beer's law[END_REF], the retrieved sources are shown in Figure5(c). Despite the reduced number of samples (41 points), the ICA technique was able to provide signals that are good approximations of the original sources. A more detailed description of this experiment is provided in[START_REF] Duarte | Application of blind source separation methods to ion-selective electrode arrays in flow-injection analysis[END_REF].

  Sources: activities of sodium and potassium.

Fig. 5 .

 5 Fig. 5. Experiment with data acquired by an array composed of a sodium and a potassium ISE.

  T denotes the source signals. The goal in source separation is to provide a good estimation of s(n). In a blind (or unsupervised) context, one has access neither to training points, i.e. a set of pairs (x i (n), s i (n)), nor to the mixing system A(•). Therefore, based only on observations x(n),

solving Equation

[START_REF] Duval | An overview of signal processing issues in chemical sensing[END_REF]

for s(n) becomes an ill-posed inverse problem whose handling requires further information (

[START_REF] Comon | Handbook of blind source separation: independent component analysis and applications[END_REF]

, chapters 7, 10, 11 and 13). For instance, a possible approach is to exploit measurement diversities such as time-frequency representations of nonstationary signals (see for example

[START_REF] Comon | Handbook of blind source separation: independent component analysis and applications[END_REF]

, chapter 11).

This work has been partially funded by the European project ERC-2012-AdG-320684-CHESS. L. T. Duarte would like to thank FAPESP (São Paulo, Brazil) for funding his research. December 19, 2013 DRAFT

However, adding this constraint can alter the separation performance if the number of sources is not correctly selected or if the mixture data are subject to attenuation and variability. A discussion on this point can be found in[START_REF] Chouzenoux | Fast constrained least squares spectral unmixing using primaldual interior-point optimization[END_REF].December 19,

DRAFT

December 19, 2013 DRAFT


Jessica Valle-Orero

Jean-Luc Garden

Jacques Richard

Andrew Wildes

Michel Peyrard

Calorimetric study of melted DNA glass

Keywords: Enthalpy Relaxation, Structural Relaxation, Deoxyribonucleic Acid, Amorphous Material, Polyvinyl Acetate, Differential Scanning Calorimetry PACS: 65.60.+a

come

INTRODUCTION

The loss of the double helical structure of DNA as a consequence of thermal denaturation has a strong effect of the whole structure of the fibres and films. While DNA is fairly rigid in its double helix form, when the two single strands dissociate after denaturation they become very flexible, and at high temperature this leads to a sample in which the strands are fully disordered and entangled, in a quasi-liquid state.

A glass transition is a transformation in which an amorphous system on cooling undergoes a transformation from liquid (equilibrium) to glass (non-equilibrium) without being able to crystalise. As the temperature of an homogeneous vitreous system in a liquid state is lowered their inter and intramolecular motions become too slow to be detectable experimentally. The liquid is impeded from reaching equilibrium which can only be attained by structural relaxation or ageing at very long timescales.

For DNA films, even though the stacking of the base pairs is broken, the intermolecular confinement is still present and influences the relaxations in the glassy state. During ageing or slow cooling the system evolves towards states with lower and lower enthalpy. Conversely, when the system is heated up after ageing, it can absorb the corresponding enthalpy at high temperature, which shows up as a peak in the heat capacity. The rate of temperature variation defines the timescales at which the system is allowed to change. The dependence of the glass structure on the cooling rate can be described by the lim-iting fictive temperature, T ′ f [1]. The glass character of DNA has been studied previously. Tsereteli's group [2] presented many published reports on the glass character of DNA in solution. Observations of the glassy behaviour in fibres were made [3] but the relaxation character of these glass-like systems was never quantitatively analysed.

In this paper, we report a detailed calorimetric study of the glassy behaviour of melted films made from DNA fibres. We compare glassy DNA melt properties with PVAc, a classical polymeric glass former. The same experimental protocols (different cooling and ageing times) are used to assess the glassy properties of the two systems. For ageing experiments, we use two different approaches concerning the evolution of the relaxation times, which lead to the same conclusion. One is based on the Tool-Narayanaswamy-Moynihan (TNM) [1] expression of the relaxation time, while the other uses an Arrhenius expression. For both approaches, the concept of limiting fictive temperature is used. It is however not usual to calculate limiting fictive temperatures after ageing experiments, it is why we focalise on this point in this paper.

EXPERIMENTAL SECTION

Sample Preparation

The samples consisted of Oriented Fibres of DNA prepared using the "wet spinning technique" which has been meticulously described elsewhere [4]. Once the oriented fibres are prepared they go through a process of cleaning, drying and humidification. For the purpose of our experiments the films were humidified to 56% RH for several weeks.

A Setaram Micro DSC III was used to perform the measurements. For each of the DSC experiments one or several hydrated films with a mass between 80 and 100 mg were rolled up and sealed inside an hastelloy tube. The sample was then placed in the DSC at 293 K and left for 10 minutes for stabilisation. For the reference tube, 345 mg of VESPEL was used, as it does not exhibit a glass transition within the region of our study.

Experimental Protocol

During the first raise in temperature from 268 K to 368 K at a rate of 1.2 K/min, the sample goes through a denaturing transition as it is observed in DNA solution. It was then cooled down to 268 K at the same rate. The cooling scan showed a glass transition at around 298 K. The sample was kept inside the DSC throughout the different measurements. Before each experiment, the sample was exposed to a series of cooling-heating temperature cycles and finally heated to 363 K during 20 mins, this is referred as "pre-experiment". This was done to erase any annealing thermal history. The data presented in specific heat △C as a function of temperature were first corrected from instrumental resolution and then the VESPEL contribution was isolated in order to obtain only the heat capacity of DNA. The resulting data and the protocol have been presented in detail in a recent publication [5].

Two series of experiments were carried out with the same sample to study the structural relaxation behaviour of the glassy system obtained by heating a film of oriented DNA fibers up to denaturation temperature and cooling it down to room temperature.

RESULTS

Cooling at different temperature rates

An enthalpic response was obtained for each of the DSC heating ramps at 1.2 K/min after cooling at the different rates β . The data expressed as normalised heat capacity as a function of temperature was published elsewhere [5]. For each of the normalised curves the limiting fictive temperature T ′ f was calculated using Equation 1 given by Moynihan [1]. The same experimental procedure and analysis was carried out for PVAc.

∫ T l T * (Cp -Cp g )dT = ∫ T l T ′ f (Cp l -Cp g )dT, ( 1
)
where Cp is the specific heat, Cp g and Cp l are the specific heat for the glassy and liquid state respectively, and T l and T * are arbitrary temperatures in the liquid and glass region further away from the transition region. The Arrhenius law,

τ ∝ 1 β ∝ exp ( △h T ′ f ) , (2)
was applied to the data in order to obtain the activation energy △h in units of temperature (in K :△h = △H * /R, where △H * is the activation enthalpy in J/mol, and R is the ideal gas constant in J/(mol•K) for both systems. These gave activation energies of 79 kK and 99 kK for DNA glass and PVAc respectively. These values are of the same order of magnitude. See Table 1.

Annealing at different ageing times

The samples were annealed at a temperature T a which is 15 K below the inflection point of their respective specific heat curve measured during the cooling run at 1.2 K/min, i.e., T a,DNA = 293 K and T a,PVAc = 297 K for DNA and PVAc. As stated before, we investigated a series of ageing times (relaxing times) t a varying from 84 h to 1 h, and even no ageing at all. In our previous work we proved that both samples behave in a qualitatively similar way [5]. However, for PVAc a peak is visible during reheating even without ageing, while for DNA there is no peak in the absence of ageing. This suggests that underlying relaxation mechanisms for the two glassformers do not obey the same laws.

The DSC heating scans for melted DNA annealed with different annealing times t a , as described in the experimental method, showed a progressive growth of an endothermic peak with the increase of t a [5]. The first step in understanding the state of the sample after annealing is to calculate the area under the endothermic peak, also called "enthalpy recovery" △H. This property corresponds to an enthalpy decrease during the annealing process, which is again absorbed by the system to reach the equilibrium state at high temperature. It characterises the annealed glassy state from the energetic viewpoint. Furthermore, the annealed system can also be characterised by another property, its limiting fictive temperature, the same way we used it in the previous section to characterise glassy states obtained at different cooling rates. Owing to our experimental protocol in which we freeze the system to low temperature (268 K) before recording the DSC signal with increasing temperature, we apply Equation 1 to the specific heat curves. This provides a value of T ′ f , which describes the instantaneous state obtained after ageing that has been frozen at lower temperature.

Once this was calculated, the data were analysed using two different expressions of the relaxation time. In both cases, the limiting fictive temperature was used for the characterisation of ageing experiments. The first method consists in a treatment based on the TNM expression of the relaxation time, in which a type of ageing time shift factor R(T a ,t a ) is used to compare the evolution during ageing of the typical relaxation times t a , between DNA and PVAc, at the ageing temperature T a [6]. The second treatment considers an Arrhenius law for the expression of the relaxation time, in which contrary to the former method, the ageing time shift factor is defined differently in order to follow the evolution of the relaxation times on the equilibrium liquid line for the two systems. Both treatments will be explained further in the following sections.

Figure 1 gives a graphical description of the concept of fictive temperature T f as a function of temperature. The graph shows how after annealing, the limiting fictive temperatures can be obtained. The continuous red curve corresponds to a cooling ramp at 1.2 K/min, where the ending point of the cooling measurement is T f = T ′ f (t a = 0), when no ageing is applied. The two characterisations of the ageing described above have been represented, the differences lie within the ending point. In the first method, the system is cooled down at the same rate, annealed at T a (marked by the arrow), and finally cooled down to a temperature lower than T a ; the relaxation time is therefore calculated at the ageing temperature τ(T a ,t a ) where the fictive temperature at that position T f can be assumed to equal T ′ f (T a ,t a ). In the second treatment, the experimental method is the same but the relaxation time is considered to be a function of the limiting fictive temperature, where T = T ′ f (T a ,t a ).

First treatment: Ratio of the relaxation times at the ageing temperature

In order to characterise the behaviour of a glass former after ageing experiments, it is usual to represent the ageing time shift factor a t a = τ(T a ,t a )/τ(T a ,t a,r ) as a function of the ageing time [6]. Where τ(T a ,t a ) is the relaxation time at the ageing temperature T a and at the ageing time t a , and t a,r is a reference ageing time. Here we present the same type of ratio using the TNM expression of the typical relaxation time:

τ = τ 0 exp ( x △h T + (1 -x) △h T f ) . (3)
In this expression, τ 0 is related to the typical relaxation time of molecular vibrations, △h is the activation enthalpy obtained by means of the different cooling rates experiments, T and T f are the temperature and the fictive temperature respectively, and x is the non-linearity parameter pertaining to the glass structure dependence of the typical relaxation time. In this case, the ageing time shift factor could be written as:

R(T a ,t a ) = τ 0 exp ( x △h T a + (1 -x) △h T f ) τ 0 exp ( △h T ′ f ) . (4)
The fictive temperature appearing in the numerator is a function of the ageing temperature and the ageing time T f = T f (T a ,t a ). From here onwards, we will assume that this fictive temperature is the same that the limiting fictive temperature calculated using Equation 1. This assumption is valid because just after the ageing time the system is cooled again to the lowest temperature before reheating. Therefore, we suppose that during this second cooling the fictive temperature of the system had no time to relax and remained in the state reached after ageing at T a . This assumption becomes more accurate as the ageing time increases. Indeed, for long ageing time the real fictive temperature of the system at T a does not have time to evolve significantly within the time interval of the second cooling. The reference relaxation time in the denominator was defined by using the limiting fictive temperature T ′ f calculated for the DSC heating ramp obtained after cooling at 1, 2 K/min. In other words, the reference ageing time is chosen for an experiment without ageing or for a reference ageing time t a,r equal to zero. The non-linearity parameter x is calculated by means of the Kovacs-Aklonis-Hutchinson-Ramos (KAHR) approach [7]. From each of the DSC heating curves obtained at the different ageing times, the enthalpy recovery (△H) and the temperature of the maximum of the peak T p were determined. Figure 2 shows T p as a function of △H for each ageing time.

The outside-saturation range shows a straight line with a rate dT p /d(△H) of about 4.03 gK/J and 1.8 gK/J for DNA and PVAc respectively. From the measurements carried out at the reference cooling rate β = 1.2 K/min, the step-like descent exhibited by the specific heat at the glass temperature △C p (T = T g ) was estimated to 0.505 J/gK for DNA and 0.52 J/gK for PVAc. Henceforth, using KHAR's method the parameter x was calculated to 0.3 for DNA and 0.48 for PVAc. Table 1 gives a summary of the characteristic values defined using the first method for both DNA and PVAc.

The data for DNA glass and PVAc was analysed using the first treatment for the ratio of the relaxation times at the ageing temperature t a . Figure 3 shows the ageing time shift factor determined as a function of the ageing temperature.

While the ageing time shift factor seems to behave similarly for DNA and PVAc for the points comprised within the 33 first hours of measurements (with only a difference in amplitude), a clear change occurs in the variation rate for the longest ageing times. The factor seems to undergo a saturation for PVAc while for DNA it increases dramatically.

Second treatment: Ratio of the relaxation times at different limiting fictive temperatures A different analysis was carried out with the same data for both systems [5]. This differs from the one described above as it is based on a different expression for the relaxation time. It aims, this time, to study the ageing dependence of the relaxation times at T ′ f , and thus to define the ageing time shift factor on the basis of the Arrhenius law. Indeed, by definition, the non-equilibrium state of a frozen system at a given temperature (in our case T a ) can be described by the corresponding equilibrium state (liquid), but taken at the temperature defined by the limiting fictive temperature. Making T f = T ′ f and T = T ′ f in the TNM expression of the relaxation time leads to an Arrhenius law. As already said, the state of the sample after annealing can also be characterised by the limiting fictive temperature. This is also consistent with the previous assumption that the fictive temperature after ageing can be frozen by the second cooling at the lowest temperature. Within this scenario, a type of ageing time shift factor can be defined such as follows: R ′ (T a ,t a ) = τ 0 exp

( △h T ′ f (T a ,t a )
) τ 0 exp

( △h T ′ f (t a =0) ) = exp ( △h T ′ f (T a ,t a ) -T ′ f (t a = 0)
) .

(5) Where T ′ f (t a = 0) is the limiting fictive temperature calculated for the reference cooling experiment with a temperature rate of 1,2 K/min (the same than for the previous treatment).

Figure 4 shows the ageing time shift factor R ′ determined as a function of the limiting fictive temperature.

In this treatment, one has to bear in mind that the evolution of the relaxation time is characteristic of the corresponding equilibrium liquid taken at different decreasing temperatures (limiting fictive temperatures), the temperature limit being the real temperature of ageing of the system. The major advantage of such a treatment is that it does not involve a complex phenomenological KAHR model.

As one can notice the ratios R and R ′ are different, but this is not surprising since they measure different quantities. R takes into account the KHAR's parameter x, while R ′ does not. And thus, the differences are not due to experimental errors. It is difficult to give a physical meaning to the absolute values of these numbers but what The straight lines corresponds to linear fits of the measurements outside the saturation range, which slope was estimated to 4.03 gK/J and 1.8 gK/J for DNA and PVAc respectively it is important is their evolution with ageing time. Both figures show that for PVAc R and R ′ vary almost linearly with ln(t a ), while for DNA glass there is a significant raise of the two ratios for long ageing times. We come back to this point in the discussion section.

DISCUSSION

Our calorimetric study of a DNA film has shown that, after the thermal denaturation of the double helix, this system has properties which are very similar to those of a polymer glass, such as PVAc, with however noticeable differences 1. The activation energies for DNA (80kK) and PVAc (100kK) are slightly different, but are on the same order of magnitude. 2. While DNA glass and PVAc seem to behave similarly for short ageing times at T g -15 K, they behave differently for ageing times larger than about 30 hours. The increase of the relaxation time for DNA seems to indicate that ageing occurs in two stages. In a first stage the system can relax rapidly to finally reach a state in which relaxation becomes very slow in the second stage. This second stage only concerns very long ageing times.

The activation energies △h for PVAc and DNA glass are measuring the energies needed to break the chemical bonds to allow the motion of the systems. For both, PVAc and DNA glass, these are weak bonds between organic elements, and thus from an energetic point of view these two compounds are not very different. This is the reason why we measure comparable activation energies (point 1). On the contrary, in point 2, we show that DNA glass relaxes slower at longer ageing times than PVAc. This difference in the relaxation times could be related to kinetic effects. For DNA, the sliding of the chains is hindered by their particular structure: they have bases dangling on the side, so that, due to the confinement of the environment, the large amplitude motions of the structure as a whole will certainly be a very slow process.

From a methodology point of view, we have shown that it is relevant to describe the ageing character of a system by the limiting fictive temperature T ′ f and hence quantify the relaxation time as a function of this property τ(T ′ f ). Similar experiments at different ageing temperatures in the glassy region could complete this study. This could bring insights into the time evolution of the typical relaxation times as a function of temperature for PVAc and DNA glass.

FIGURE 1 .

 1 FIGURE 1. Evolution of the fictive temperature versus temperature after annealing; see text for details.

FIGURE 2 .

 2 FIGURE 2. Dependence of the recovery enthalpy △H on the temperature of the enthalpic peak T p for DNA (a) and PVAc (b).The straight lines corresponds to linear fits of the measurements outside the saturation range, which slope was estimated to 4.03 gK/J and 1.8 gK/J for DNA and PVAc respectively

FIGURE 3 .

 3 FIGURE 3. Evolution of the ageing time shift factor R(T a ,t a ) versus the logarithm of the ageing time ln(t a ) for DNA (open circles) and PVAc (filled circles). The data was analysed using the first treatment, in which the relaxation time is defined at the ageing temperature.

FIGURE 4 .

 4 FIGURE 4. Evolution of the ageing time shift factor R ′ (T a ,t a ) versus the logarithm of the ageing time ln(t a ) for DNA (open circles) and PVAc (filled circles). The data was analysed using the second treatment, in which the relaxation time is defined at the limiting fictive temperature T ′ f .

  Annealing at different ageing times. The procedure consisted in cooling at 1.2 K/min from 363 K to the ageing temperature T a below the glass transition temperature, annealing the sample for a period of time t a , and then cooling down to 268 K before recording the DSC signal on heating from 268 K to 363 K. All temperature ramps were programmed at a rate of 1.2 K/min. The relaxation of the system was studied at fifteen different ageing times, t a = 84, 78, 72, 66, 60, 54, 48, 42, 32, 20, 10, 6, 3, 1 and 0 hours. The last corresponds to an experiment without annealing. However by simply taking into account the minimal response time of the calorimeter one can consider that the effective minimal ageing time is actually equal to 10 min.

	ing rate was kept fixed.
	2.
	1. Different cooling rates. This consisted of a succes-
	sion of temperature ramps. Starting from the initial
	state at 363 K, the sample was cooled down at a rate
	β from 363 to 268 K; left for 5 mins for equilibra-
	tion; and heated to 363 K at 1.2 K/min. It was then
	stabilised for 20 min at this highest temperature be-
	fore the next cooling ramp, which was recorded as
	a DSC scan. The cooling rate was modified at each
	ramp in the following order, β = 1.2, 0.4, 0.2, 0.1,
	0.05, 0.033, 0.01 and back to 1.2 K/min. The heat-

TABLE 1 .

 1 Summary of the calculated values relevant for the analysis of the relaxation behaviour of both systems, DNA and PVAc. The first two, △H and △C p , were necessary in order to obtain the non-linearity parameter x.T a (K) T g (K)

					dT p
				dT p d△H (gK/J) △C p (J/gK)	x	△h (kK) Molecular weight
	DNA	293	308	4.03	0.51	0.30	79	1.3 x10 7
	PVAc	297	312	1.80	0.52	0.48	99	1.6 x 10 5