Andrej Pázman
email: pazman@fmph.uniba.sk

Luc Pronzato
email: pronzato@i3s.unice.fr

OPTIMUM DESIGN ACCOUNTING FOR THE GLOBAL NONLINEAR BEHAVIOR OF THE MODEL

Keywords: optimal design, nonlinear least-squares, estimability, curvature

Among the major difficulties that one may encounter when estimating parameters in a nonlinear regression model are the nonuniqueness of the estimator, its instability with respect to small perturbations of the observations and the presence of local optimizers of the estimation criterion.

We show that these estimability issues can be taken into account at the design stage, through the definition of suitable design criteria.

Extensions of E-, c-and G-optimality criteria are considered, which, when evaluated at a given θ 0 (local optimal design), account for the behavior of the model response η(θ) for θ far from θ 0 . In particular, they ensure some protection against close-to-overlapping situations where η(θ) -η(θ 0 ) is small for some θ far from θ 0 . These extended criteria are concave and necessary and sufficient conditions for optimality (Equivalence Theorems) can be formulated. They are not differentiable, but when the design space is finite and the set Θ of admissible θ is discretized, optimal design forms a linear programming problem which can be solved directly or via relaxation when Θ is just compact. Several examples are presented.

1. Introduction. We consider a nonlinear regression model with observations y i = y(x i ) = η(x i , θ) + ε i , i = 1, . . . , N ,

where the errors ε i satisfy IE(ε i ) = 0, var(ε i ) = σ 2 and cov(ε i , ε j ) = 0 for i = j, i, j = 1, . . . , N , and the true value θ of the vector of model parameter θ belongs to Θ, a compact subset of R p such that Θ ⊂ int(Θ), the closure of the interior of Θ. In a vector notation, we write

(1)

y = η X ( θ) + ε , with IE(ε) = 0 , Var(ε) = σ 2 I N ,
where η X (θ) = (η(x 1 , θ), . . . , η(x N , θ)) ⊤ , y = (y 1 , . . . , y N ) ⊤ , ε = (ε 1 , . . . , ε N ) ⊤ , and X denotes the N -point exact design (x 1 , . . . , x N ). The more general nonstationary (heteroscedastic) case where var(ε i ) = σ 2 (x i ) can easily be transformed into the model (1) with σ 2 = 1 via the division of y i and η(x i , θ) by σ(x i ). We suppose that η(x, θ) is twice continuously differentiable with respect to θ ∈ int(Θ) for any x ∈ X , a compact subset of R d . The model is assumed to be identifiable over X ; that is, we suppose that

(2) η(x, θ ′ ) = η(x, θ) for all x ∈ X =⇒ θ ′ = θ .

We shall denote by Ξ the set of design measures ξ, i.e., of probability measures on X . The information matrix (for σ 2 = 1) for the design X at θ is

M(X, θ) = N i=1 ∂η(x i , θ) ∂θ ∂η(x i , θ) ∂θ ⊤
and, for any ξ ∈ Ξ, we shall write

M(ξ, θ) = X [∂η(x, θ)/∂θ] [∂η(x, θ)/∂θ ⊤ ] ξ(dx) .
Denoting ξ N = (1/N ) N i=1 δ x i the empirical design measure associated with X, with δ x the delta measure at x, we have M(X, θ) = N M(ξ N , θ). Note that (2) implies the existence of a ξ ∈ Ξ satisfying the Least-Squares (LS) estimability condition

(3) η(x, θ ′ ) = η(x, θ) ξ-almost everywhere =⇒ θ ′ = θ .

Given an exact N -point design X, the set of all hypothetical means of the observed vectors y in the sample space R N forms the expectation surface S η = {η X (θ) : θ ∈ Θ}. Since η X (θ) is supposed to have continuous first and second-order derivatives in int(Θ), S η is a smooth surface in R N with a (local) dimension given by r = rank[∂η X (θ)/∂θ ⊤ ]. If r = p (which means full rank), the model ( 1) is said regular. In regular models with no overlapping of S η , i.e. when η X (θ) = η X (θ ′ ) implies θ = θ ′ , the LS estimator (4) θLS = θN LS = arg min θ∈Θ y -η X (θ) 2 is uniquely defined with probability one (w.p.1). Indeed, when the distributions of errors ε i have probability densities (in the standard sense) it can be proven that η[ θLS (y)] is unique w.p.1, see [START_REF] Pázman | Nonlinear least squares -Uniqueness versus ambiguity[END_REF] and Pázman (1993, p. 107). However, there is still a positive probability that the function θ -→ y -η X (θ) 2 has a local minimizer different from the global one when the regression model is intrinsically curved in the sense of [START_REF] Bates | Relative curvature measures of nonlinearity[END_REF], i.e., when S η is a curved surface in R N , see [START_REF] Demidenko | Optimization and Regression[END_REF][START_REF] Demidenko | Is this the least squares estimate?[END_REF]. Moreover, a curved surface can "almost overlap"; that is, there may exist points θ and θ ′ in Θ such that θ ′ -θ is large but η X (θ ′ ) -η X (θ) is small (or even equals zero in case of strict overlapping). This phenomenon can cause serious difficulties in parameter estimation, leading to instabilities of the estimator, and one should thus attempt to reduce its effects by choosing an adequate experimental design. Classically, those issues are ignored at the design stage and the experiment is chosen on the basis of asymptotic local properties of the estimator. Even when the design relies on small-sample properties of the estimator, like in [START_REF] Pázman | Nonlinear experimental design based on the distribution of estimators[END_REF][START_REF] Gauchi | Designs in nonlinear regression by stochastic minimization of functionnals of the mean square error matrix[END_REF], a non-overlapping assumption is used (see Pázman (1993, pp. 66 and 157)) which permits to avoid the aforementioned difficulties. Note that putting restrictions on curvature measures is not enough: consider the case dim(θ) = 1 with the overlapping S η formed by a circle of arbitrarily large radius and thus arbitrarily small curvature (see the example in Sect. 2 below).

Important and precise results are available concerning the construction of subsets of Θ where such difficulties are guaranteed not to occur, see, e.g., [START_REF] Chavent | Local stability of the output least square parameter estimation technique[END_REF][START_REF] Chavent | A new Sufficient Condition for the Wellposedness of Non-Linear Least-Square Problems arising in Identification and Control[END_REF][START_REF] Chavent | New Size × Curvature Conditions for Strict Quasiconvexity of Sets[END_REF]; however, their exploitation for choosing adequate designs is far from straightforward. Also, the construction of designs with restricted curvatures, as proposed by [START_REF] Chaloner | Constrained design strategies for improving normal approximations in nonlinear regression problems[END_REF], is based on the curvature measures of [START_REF] Bates | Relative curvature measures of nonlinearity[END_REF] and uses derivatives of η X (θ) at a certain θ; this local approach is unable to catch the problem of overlapping for two points that are distant in the parameter space. Other design criteria using a second-order development of the model response, or an approximation of the density of θLS [START_REF] Hamilton | A quadratic design criterion for precise estimation in nonlinear regression models[END_REF][START_REF] Pronzato | Second-order approximation of the entropy in nonlinear least-squares estimation[END_REF], are also inadequate.

The aim of this paper is to present new optimality criteria for optimum design in nonlinear regression models that may reduce such effects, especially overlapping, and are at the same time closely related to classical optimality criteria like E-, c-or G-optimality (in fact, they coincide with those criteria when the regression model is linear). Classical optimality criteria focus on efficiency, i.e., aim at ensuring a precise estimation of θ, asymptotically, provided that the model is locally identifiable at θ. On the other hand, the new extended criteria account for the global behavior of the model and enforce identifiability.

An elementary example is given in the next section and illustrates the motivation of our work. The criterion of extended E-optimality is considered in Sect. 3; its main properties are detailed and algorithms for the construction of optimal designs are presented. Sections 4 and 5 are respectively devoted to the criteria of extended c-optimality and extended G-optimality. Several illustrative examples are presented in Sect. 6. Section 7 suggests some extensions and further developments and Sect. 8 concludes.

2. An elementary motivating example.

Example 1. Suppose that θ ∈ Θ = [0, 1] and that, for any design point

x = (t, u) ⊤ ∈ X = {0, π/2} × [0, u max ], we have η(x, θ) = r cos(t -u θ) ,
with r a known positive constant. We take u max = 7 π/4; the difficulties mentioned below are even more pronounced for values of u max closer to 2π. We shall consider exclusively two-point designs X = (x 1 , x 2 ) of the form

x 1 = (0, u) ⊤ , x 2 = (π/2, u) ⊤
and denote ν u the associated design measure,

ν u = (1/2)[δ x 1 + δ x 2 ].
We shall look for an optimal design, that is, an optimal choice of u ∈ [0, u max ], where optimality is considered in terms of information.

It is easy to see that for any design ν u we have

η X (θ) = η(x 1 , θ) η(x 2 , θ) = r cos(u θ) r sin(u θ) .
The expectation surface is then an arc of a circle, with central angle u, see Figure 1 for the case u = u max = 7π/4. The model is nonlinear but parametrically linear since the information matrix M (X, θ) for σ 2 = 1 (here scalar since θ is scalar) equals r 2 u 2 and does not depend on θ. Also, the intrinsic curvature (see Sect. 6) is constant and equals 1/r, and the model is also almost intrinsically linear if r gets large. Any classical optimality criterion (A-, D-, E-) indicates that one should observe at u = u max , and setting a constraint on the intrinsic curvature is not possible here. However, if the true value of θ is θ = 0 and σ 2 is large enough, there is a chance that the LS estimator will be θLS = 1, and thus very far from θ, see Figure 1. The situation gets even worse if u max gets closer to 2 π, since S η then almost overlaps. Now, consider expectation surface S η is then a half-circle, so that η X (0) and η X (1) are as far away as possible. On the other hand, as shown in Sect. 3, min θ∈Θ H E (ν u , θ) possesses most of the attractive properties of classical optimality criteria and even coincides with one of them in linear models. Figure 2-left shows H E (ν u , θ) as a function of θ for three values of u and illustrates the fact that the minimum of H E (ν u , θ) with respect to θ ∈ Θ is maximized for u = u * . Figure 2-right shows that the design with u = u max (dashed line) is optimal locally at θ = θ 0 , in the sense that it yields the fastest increase of η X (θ) -η X (θ 0 ) as θ slightly deviates from θ 0 . On the other hand, u = π maximizes min θ∈Θ η X (θ) -η X (θ 0 ) /|θ -θ 0 | (solid line) and realizes a better protection against the folding effect of S η , at the price of a slightly less informative experiment for θ close to θ 0 . Smaller values of u (dotted line) are worse than u * , both locally for θ close to θ 0 and globally in terms of the folding of S η .

H E (ν u , θ) = (1/2) η X (θ)-η X (θ 0 ) 2 /|θ -θ 0 | 2 , see (5), with θ 0 = 0. For all u ∈ [0, u max ], the minimum of H E (ν u , θ) with respect to θ ∈ Θ is obtained at θ = 1, H E (ν u , 1) = r 2 [1-cos(u)] is then maximum in [0, u max ] for u = u * = π.
The rest of the paper will formalize these ideas and show how to implement them for general nonlinear models through the definition of suitable design criteria that can be easily optimized.

Extended (globalized) E-optimality.

3.1. Definition of φ eE (•). Take a fixed point θ 0 in Θ and denote (5) where • ξ denotes the norm in L 2 (ξ); that is, l ξ = X l 2 (x) ξ(dx)

H E (ξ, θ) = H E (ξ, θ; θ 0 ) = η(•, θ) -η(•, θ 0 ) 2 ξ θ -θ 0 2 , 0 0.
1/2 for any l ∈ L 2 (ξ). When ξ is a discrete measure, like in the examples considered in the paper, then l 2 ξ is simply the sum x:ξ({x})>0 ξ({x}) l 2 (x). The extended E-optimality criterion is defined by

(6) φ eE (ξ) = φ eE (ξ; θ 0 ) = min θ∈Θ H E (ξ, θ) ,
to be maximized with respect to the design measure ξ.

In a nonlinear regression model φ eE (•) depends on the value chosen for θ 0 and can thus be considered as a local optimality criterion. On the other hand, the criterion is global in the sense that it depends on the behavior of η(•, θ) for θ far from θ 0 . This (limited) locality can be removed by considering φ M eE (ξ) = min θ 0 ∈Θ φ eE (ξ; θ 0 ) instead of ( 6), but only the case of φ eE (•) will be detailed in the paper, the developments being similar for φ M eE (•), see Sect. 7.2.

For a linear regression model with η(x, θ) = f ⊤ (x)θ+v(x) and Θ = R p , for any θ 0 and any ξ ∈ Ξ we have η(

•, θ) -η(•, θ 0 ) 2 ξ = (θ -θ 0 ) ⊤ M(ξ)(θ -θ 0 ), so that φ eE (ξ) = min θ-θ 0 ∈R p (θ -θ 0 ) ⊤ M(ξ)(θ -θ 0 ) θ -θ 0 2 = λ min [M(ξ)] ,
the minimum eigenvalue of M(ξ), and corresponds to the E-optimality criterion.

For a nonlinear model with Θ = B(θ 0 , ρ), the ball with centre θ 0 and radius ρ, direct calculation shows that In a nonlinear regression model with larger Θ, the determination of an optimum design ξ * eE maximizing φ eE (ξ) ensures some protection against η(•, θ) -η(•, θ 0 ) ξ being small for some θ far from θ 0 . In particular, when θ 0 ∈ int(Θ) then φ eE (ξ; θ 0 ) = 0 if either M(ξ, θ 0 ) is singular or η(•, θ)η(•, θ 0 ) ξ = 0 for some θ = θ 0 . Therefore, under the condition (2), ξ * eE satisfies the estimability condition (3) at θ = θ 0 and is necessarily nondegenerate, i.e., M(ξ * eE , θ 0 ) is nonsingular, when θ 0 ∈ int(Θ) (provided that there exists a non-degenerate design in Ξ). Notice that (7) implies that φ eE (ξ; θ 0 ) ≤ λ min [M(ξ, θ 0 )] when Θ contains some open neighborhood of θ 0 . In contrast with the E-optimality criterion, maximizing φ eE (ξ; θ 0 ) in nonlinear models does not require computation of the derivatives of η(x, θ) with respect to θ at θ 0 , see the algorithms proposed in Sects. 3.3 and 3.4. Also note that the influence of points that are very far from θ 0 can be suppressed by modification of the denominator of (5) without changing the relation with E-optimality, see Sect. 7.1.

Before investigating properties of φ eE (•) as a criterion function for optimum design in the next section, we state a property relating φ eE (ξ) to the localization of the LS estimator θLS .

Theorem 1 For any given θ ∈ Θ, the LS estimator θLS given by (4) in the model (1) satisfies

θLS ∈ Θ ∩ B θ, 2 y -η X (θ) √ N φ eE (ξ N ; θ) ,
with ξ N the empirical measure associated with the design X used to observe y.

Proof. The result follows from the following chain of inequalities

θLS -θ ≤ η(•, θLS ) -η(•, θ) ξ N φ eE (ξ N ; θ) = η X ( θLS ) -η X (θ) √ N φ eE (ξ N ; θ) ≤ y -η X ( θLS ) + y -η X (θ) √ N φ eE (ξ N ; θ) ≤ 2 y -η X (θ) √ N φ eE (ξ N ; θ) . (8)
Note that although the bound ( 8) is tight in general nonlinear situations (due to the possibility that S η overlaps), it is often pessimistic. In particular, in the linear regression model η

(x, θ) = f ⊤ (x)θ+v(x), direct calculation gives θLS -θ ≤ λ max [(F ⊤ F) -1 ] y -η X (θ) = y -η X (θ) √ N φ eE (ξ N ) ,
where F is the N × p matrix with ith line equal to f ⊤ (x i ). We also have θLS

-θ ≤ y -η X (θ) /[ √ N φ eE (ξ N , θ)] in intrinsically linear models (with a flat expectation surface S η ) since then η X ( θLS ) -η X (θ) ≤ y - η X (θ) .
In the following we shall omit the dependence in θ 0 and simply write φ eE (ξ) for φ eE (ξ; θ 0 ) when there is no ambiguity.

Properties of φ eE (•). As the minimum of linear functions of

ξ, φ eE (•) is concave: for all ξ, ν ∈ Ξ and all α ∈ [0, 1], φ eE [(1 -α)ξ + αν] ≥ (1 - α)φ eE (ξ) + αφ eE (ν). It is also positively homogeneous: φ eE (aξ) = a φ eE (ξ)
for all ξ ∈ Ξ and a > 0, see, e.g., Pukelsheim (1993, Chap. 5). The criterion of eE-efficiency can then be defined as

E eE (ξ) = φ eE (ξ) φ eE (ξ * eE ) , ξ ∈ Ξ ,
where ξ * eE maximizes φ eE (ξ). The concavity of φ eE (•) implies the existence of directional derivatives and, due to the linearity in ξ of H E (ξ, θ), we have the following, see, e.g., [START_REF] Dem'yanov | Introduction to Minimax[END_REF].

Theorem 2 For any ξ, ν ∈ Ξ, the directional derivative of the criterion φ eE (•) at ξ in the direction ν is given by

F φ eE (ξ; ν) = min θ∈Θ E (ξ) H E (ν, θ) -φ eE (ξ) ,
where

Θ E (ξ) = {θ ∈ Θ : H E (ξ, θ) = φ eE (ξ)}. Note that we can write F φ eE (ξ; ν) = min θ∈Θ E (ξ) X Ψ eE (x, θ, ξ) ν(dx), where (9) Ψ eE (x, θ, ξ) = [η(x, θ) -η(x, θ 0 )] 2 -η(•, θ) -η(•, θ 0 ) 2 ξ θ -θ 0 2 .
Due to the concavity of φ eE (•), a necessary and sufficient condition for the optimality of a design measure ξ * eE is that

(10) sup ν∈Ξ F φ eE (ξ * eE ; ν) = 0 ,
a condition often called "Equivalence Theorem" in optimal design theory; see, e.g., [START_REF] Fedorov | Theory of Optimal Experiments[END_REF]; [START_REF] Silvey | Optimal Design[END_REF]. An equivalent condition is as follows.

imsart-aos ver. 2013/03/06 file: Extended-optimality_AOS-REV.tex date: March 31, 2014

Theorem 3 A design ξ * eE ∈ Ξ is optimal for φ eE (•) if and only if (11) max x∈X Θ E (ξ * eE ) Ψ eE (x, θ, ξ) µ * (dθ) = 0 for some measure µ * ∈ M [Θ E (ξ * eE )] ,
the set of probability measures on Θ E (ξ * eE ). Proof. This is a classical result for maximin design problems, see, e.g., Fedorov and Hackl (1997, Sect. 2.6). We have

0 ≤ sup ν∈Ξ F φ eE (ξ; ν) = sup ν∈Ξ min θ∈Θ E (ξ) X Ψ eE (x, θ, ξ) ν(dx) = sup ν∈Ξ min µ∈M [Θ E (ξ)] X Θ E (ξ) Ψ eE (x, θ, ξ) µ(dθ) ν(dx) = min µ∈M [Θ E (ξ)] sup ν∈Ξ X Θ E (ξ) Ψ eE (x, θ, ξ) µ(dθ) ν(dx) = min µ∈M [Θ E (ξ)] max x∈X Θ E (ξ) Ψ eE (x, θ, ξ) µ(dθ) . (12)
Therefore, the necessary and sufficient condition (10) can be written as (11).

One should notice that sup ν∈Ξ F φ eE (ξ; ν) is generally not obtained for ν equal to a one-point (delta) measure, which prohibits the usage of classical vertex-direction algorithms for optimizing φ eE (•). Indeed, the minimax problem (12) has generally several solutions x (i) for x, i = 1, . . . , s, and the optimal ν * is then a linear combination s i=1 w i δ x (i) , with w i ≥ 0 and s i=1 w i = 1; see [START_REF] Pronzato | Nonsequential T -optimal design for model discrimination: new algorithms[END_REF] for developments on a similar difficulty in T -optimum design for model discrimination. This property, due to the fact that φ eE (•) is not differentiable, has the important consequence that the determination of a maximin-optimal design cannot be obtained via standard design algorithms used for differentiable criteria.

To avoid that difficulty, a regularized version φ eE,λ (•) of φ eE (•) is considered in (Pronzato and Pázman, 2013, Sects. 7.7.3 and 8.3.2), with the property that lim λ→∞ φ eE,λ (ξ) = φ eE (ξ) for any ξ ∈ Ξ (the convergence being uniform when Θ is a finite set), φ eE (•) is concave and such that sup ν∈Ξ F φ eE,λ (ξ; ν) is obtained when ν is the delta measure δ x * at some x * ∈ X (depending on ξ). However, although φ eE,λ (•) is smooth for any finite λ, its maximization tends to be badly conditioned for large λ.

In the next section we show that optimal design for φ eE (•) reduces to linear programming when Θ and X are finite. This is an important property. An algorithm based on a relaxation of the maximin problem is then considered in Sect. 3.4 for the case where Θ is compact.

3.3. Optimal design via linear-programming (Θ is finite). To simplify the construction of an optimal design, one may take Θ as a finite set, Θ = Θ (m) = {θ (1) , θ (2) , . . . , θ (m) }; φ eE (ξ) can then be written as φ eE (ξ) = min j=1,...,m H E (ξ, θ (j) ), with H E (ξ, θ) given by ( 5). If the design space X is also finite, with X = {x (1) , x (2) , . . . , x (ℓ) }, then the determination of an optimal design measure for φ eE (•) amounts to the determination of a scalar t and of a vector of weights w = (w 1 , w 2 , . . . , w ℓ ) ⊤ , w i being allocated at x (i) for each i = 1, . . . , ℓ, such that c ⊤ [w ⊤ , t] ⊤ is maximized, with c = (0, 0, . . . , 0, 1) ⊤ and w and t satisfying the constraints

ℓ i=1 w i = 1 , w i ≥ 0 , i = 1, . . . , ℓ , ℓ i=1 w i h i (θ (j) ) ≥ t , j = 1, . . . , m , (13)
where we denoted ( 14)

h i (θ) = [η(x (i) , θ) -η(x (i) , θ 0 )] 2 θ -θ 0 2 .
This is a linear programming (LP) problem, which can easily be solved using standard methods (the simplex algorithm for instance), even for large m and ℓ. We shall denote by ( ŵ, t) = LP eE (X , Θ (m) ) the solution of this problem.

We show below how a compact subset Θ of R p with nonempty interior can be replaced by a suitable discretized version Θ (m) that can be enlarged iteratively.

3.4. Optimal design via relaxation and the cutting-plane method (Θ is a compact subset of R p ). Suppose now that X is finite and that Θ is a compact subset of R p with nonempty interior. In the LP formulation above, (w, t) must satisfy an infinite number of constraints: ℓ i=1 w i h i (θ) ≥ t for all θ ∈ Θ, see (13). One may then use the method of [START_REF] Shimizu | Necessary conditions for min-max problems and algorithm by a relaxation procedure[END_REF] and consider the solution of a series of relaxed LP problems, using at step k a finite set of constraints only, i.e., consider θ ∈ Θ (k) finite. Once a solution (w k , t k ) = LP eE (X , Θ (k) ) of this problem is obtained, using a standard LP solver, the set Θ (k) is enlarged to Θ (k+1) = Θ (k) ∪ {θ (k+1) } with imsart-aos ver. 2013/03/06 file: Extended-optimality_AOS-REV.tex date: March 31, 2014 θ (k+1) given by the constraint (13) most violated by w k , i.e., ( 15)

θ (k+1) = arg min θ∈Θ H E (w k , θ) ,
where, with a slight abuse of notation, we write H E (w, θ) = H E (ξ, θ), see (5), when ξ allocates mass w i at the support point x (i) ∈ X for all i. This yields the following algorithm for the maximization of φ eE (•). 0) Take any vector w 0 of nonnegative weights summing to one, choose ǫ > 0, set Θ (0) = ∅ and k = 0. 1) Compute θ (k+1) given by ( 15), set Θ

(k+1) = Θ (k) ∪ {θ (k+1) }. 2) Use a LP solver to determine (w k+1 , t k+1 ) = LP eE (X , Θ (k+1) ) 3) If ∆ k+1 = t k+1 -φ eE (w k+1
) < ǫ, take w k+1 as an ǫ-optimal solution and stop; otherwise k ← k + 1, return to step 1.

The optimal value φ * eE = max ξ∈Ξ φ eE (ξ) satisfies

φ eE (w k+1 ) ≤ φ * eE ≤ t k+1
at every iteration, so that ∆ k+1 of step 3 gives an upper bound on the distance to the optimum in terms of criterion value. The algorithm can be interpreted in terms of the cutting-plane method. Indeed, from ( 5) and ( 14) we have

H E (w, θ (j+1) ) = ℓ i=1 w i h i (θ (j+1)
) for any vector of weights w. From the definition of θ (j+1) in ( 15) we obtain

φ eE (w) ≤ H E (w, θ (j+1) ) = H E (w j , θ (j+1) ) + ℓ i=1 h i (θ (j+1) ){w -w j } i = φ eE (w j ) + ℓ i=1 h i (θ (j+1) ){w -w j } i ,
so that the vector with components h i (θ (j+1) ), i = 1, . . . , ℓ, forms a subgradient of φ eE (•) at w j , which we denote ∇φ eE (w j ) below (it is sometimes called supergradient since φ eE (•) is concave). Each of the constraints

ℓ i=1 w i h i (θ (j+1) ) ≥ t ,
used in the LP problem of step 2, with j = 0, . . . , k, can be written as Therefore, w k+1 determined at step 2 maximizes the piecewise-linear approximation min j=0,...,k {φ eE (w j ) + ∇ ⊤ φ eE (w j )(ww j )} of φ eE (w) with respect to the vector of weights w, and the algorithm corresponds to the cutting-plane method of [START_REF] Kelley | The cutting plane method for solving convex programs[END_REF].

∇ ⊤ φ eE (w j )w = φ eE (w j ) + ∇ ⊤ φ eE (w j )(w -w j ) ≥ t .
The only difficult step in the algorithm corresponds to the determination of θ (k+1) in ( 15) when Θ is a compact set. We found that the following simple procedure is rather efficient. Construct a finite grid, or a space-filling design, G 0 in Θ. Then, for k = 0, 1, 2, . . .

(16)              (i) compute θk+1 = arg min θ ′ ∈G k H E (w k , θ ′ ) ; (ii) perform a local minimization of H E (w k , θ) with respect to θ ∈ Θ , initialized at θk+1 ; let θ (k+1) denote the solution ; (iii) set G k+1 = G k ∪ {θ (k+1) } .
The optimal value φ eE (ξ * eE ) can then be approximated by H E (w k+1 , θ (k+2) ) when the algorithm stops (step 3).

The method of cutting planes is known to have sometimes rather poor convergence properties, see, e.g., Bonnans et al. (2006, Chap. 9), Nesterov (2004, Sect. 3.3.2). A significant improvement consists in restricting the search for w k+1 at step 2 to some neighborhood of the best solution obtained so far, which forms the central idea of bundle methods, see [START_REF] Lemaréchal | New variants of bundle methods[END_REF], Bonnans et al. (2006, Chaps. 9-10). In particular, the level method of Nesterov (2004, Sect. 3.3.3) adds a quadratic-programming step to each iteration of the cutting planes algorithm presented above; one may refer for instance to Pronzato and Pázman (2013, Sect. 9.5.3) for an application of the level method to design problems. Notice that any linear constraint on w can easily be taken into account in addition to those in (13), so that the method directly applies to optimal design with linear cost-constraints, see, e.g., Fedorov and Leonov (2014, Sect. 4.2).

4. Extended (globalized) c-optimality.

H c (ξ, θ) = H c (ξ, θ; θ 0 ) = η(•, θ) -η(•, θ 0 ) 2 ξ |g(θ) -g(θ 0 )| 2
and consider the design criterion defined by

(18) φ ec (ξ) = min θ∈Θ H c (ξ, θ) ,
to be maximized with respect to the design measure ξ.

When η(x, θ) and the scalar function g(θ) are both linear in θ, with g(θ) = c ⊤ θ, we get Harville (1997, eq. 10.4). Also, for a nonlinear model with Θ = B(θ 0 , ρ) and a design ξ such that M(ξ, θ 0 ) has full rank, one has lim

φ ec (ξ) = min θ∈Θ, c ⊤ (θ-θ 0 ) =0 (θ -θ 0 ) ⊤ M(ξ)(θ -θ 0 ) [c ⊤ (θ -θ 0 )] 2 and therefore φ ec (ξ) = [c ⊤ M -(ξ)c] -1 , using the well-known formula c ⊤ M -c = max α =0 (c ⊤ α) 2 /(α ⊤ Mα), c.f.
ρ→0 φ ec (ξ) = [c ⊤ M -1 (ξ, θ 0 )c] -1 ,
which justifies that we consider φ ec (ξ) as an extended c-optimality criterion.

At the same time, in a nonlinear situation with larger Θ the determination of an optimal design ξ * ec maximizing φ ec (ξ) ensures some protection against η(•, θ) -η(•, θ 0 ) 2 ξ being small for some θ such that g(θ) is significantly different from g(θ 0 ). The condition (2) guarantees the existence of a ξ ∈ Ξ such that φ ec (ξ) > 0 and thus the LS estimability of g(θ) at θ 0 for ξ * ec , that is, η(x, θ) = η(x, θ 0 ) ξ * ec -almost everywhere =⇒ g(θ) = g(θ 0 ) , see Pronzato and Pázman (2013, Sect. 7.4.4). When Θ contains an open neighborhood of θ 0 , then

φ ec (ξ) ≤ [c ⊤ M -(ξ, θ 0 )c] -1 .
Similarly to φ eE (•), the criterion φ ec (•) is concave and positively homogeneous; its concavity implies the existence of directional derivatives.

Theorem 4 For any ξ, ν ∈ Ξ, the directional derivative of the criterion φ ec (•) at ξ in the direction ν is given by A necessary and sufficient condition for the optimality of ξ * maximizing φ ec (•) is that sup ν∈Ξ F φec (ξ * ; ν) = 0, which yields an Equivalence Theorem similar to Th. 3.

F φec (ξ; ν) = min θ∈Θc(ξ) H c (ν, θ) -φ ec (ξ) , where Θ c (ξ) = {θ ∈ Θ : H c (ξ, θ) = φ ec (ξ)}.
When both Θ and X are finite, an optimal design for φ ec (•) is obtained by solving a LP problem. Compared with Sect. 3.3, we simply need to substitute H c for H E and use h i (θ) = [η(x (i) , θ) -η(x (i) , θ 0 )] 2 /|g(θ) -g(θ 0 )| 2 , i = 1, . . . , ℓ, instead of ( 14). Also, a relaxation method similar to that in Sect. 3.4 can be used when Θ is a compact subset of R p .

5. Extended (globalized) G-optimality. Following the same lines as above, we can also define an extended G-optimality criterion by

φ eG (ξ) = min θ∈Θ η(•, θ) -η(•, θ 0 ) 2 ξ max x∈X [η(x, θ) -η(x, θ 0 )] 2 .
The fact that it corresponds to the G-optimality criterion for a linear model can easily be seen, noticing that in the model (1) with η(x, θ) = f ⊤ (x)θ+v(x) we have

sup x∈X N σ 2 var f ⊤ (x) θLS -1 = inf x∈X [f ⊤ (x)M -1 (ξ N )f (x)] -1 = inf x∈X inf u∈R p , u ⊤ f (x) =0 u ⊤ M(ξ N )u [f ⊤ (x)u] 2 = inf u∈R p u ⊤ M(ξ N )u max x∈X [f ⊤ (x)u] 2 ,
where ξ N denotes the empirical design measure corresponding to X, assumed to be nonsingular, and the second equality follows from Harville (1997, eq. 10.4). The equivalence theorem of [START_REF] Kiefer | The equivalence of two extremum problems[END_REF] indicates that D-and G-optimal designs coincide; therefore, D-optimal designs are optimal for φ eG (•) in linear models. Moreover, the optimum (maximum) value of φ eG (ξ) equals 1/p with p = dim(θ).

In a nonlinear model, a design ξ * eG maximizing φ eG (ξ) satisfies the estimability condition (3) at θ = θ 0 . Indeed, max x∈X η(x, θ) -η(x, θ 0 ) 2 > 0 for any θ = θ 0 from (2), so that there exists some ξ ∈ Ξ such that φ eG (ξ) > 0. Therefore, φ eG (ξ * eG ) > 0, and η(•, θ)-η(•, θ 0 ) 2 ξ * eG = 0 implies that η(x, θ) = η(x, θ 0 ) for all x ∈ X , that is, θ = θ 0 from (2). Notice that when Θ contains an open neighborhood of θ 0 , then φ eG (ξ) ≤ 1/p for all ξ ∈ Ξ.

Again, directional derivatives can easily be computed and an optimal design can be obtained by linear programming when Θ and X are both finite, or with the algorithm of Sect. 3.4 when X is finite but Θ has nonempty imsart-aos ver. 2013/03/06 file: Extended-optimality_AOS-REV.tex date: March 31, 2014 interior. Note that there are now m × ℓ inequality constraints in (13), given by ℓ i=1 w i h i (θ (j) , x (k) ) ≥ t , j = 1, . . . , m , k = 1, . . . , ℓ , where now

h i (θ, x) = [η(x (i) , θ) -η(x (i) , θ 0 )] 2
[η(x, θ) -η(x, θ 0 )] 2 . Also note that in the algorithm of Sect. 3.4 we need to construct two sequences of sets, Θ (k) and X (k) , with Θ (k+1) = Θ (k) ∪ {θ (k+1) } and X (k+1) = X (k) ∪ {x (k+1) } at step 2, and (15) replaced by

{θ (k+1) , x(k+1) } = arg min {θ,x}∈Θ×X η(•, θ) -η(•, θ 0 ) 2 ξ k [η(x, θ) -η(x, θ 0 )] 2
with ξ k the design measure corresponding to the weights w k .

6.

Examples. We shall use the common notation

ξ = x 1 • • • x m w 1 • • • w m
for a discrete design measure with m support points x i and such that ξ({x i }) = w i , i = 1, . . . , m. In the three examples considered we indicate the values of the parametric, intrinsic and total measure of curvatures at θ 0 (for σ = 1), see Tables 1, 2 and 3. They are not used for the construction of optimal designs, and the examples illustrate the fact that they provide information on the local behavior only (at θ 0 ), so that a small curvature does not mean good performance in terms of extended optimality. They are given by

C int (ξ, θ) = sup u∈R p -{0} [I -P θ ] p i,j=1 u i [∂ 2 η(•, θ)/∂θ i ∂θ j ]u j ξ u ⊤ M(ξ, θ)u , C par (ξ, θ) = sup u∈R p -{0} P θ p i,j=1 u i [∂ 2 η(•, θ)/∂θ i ∂θ j ]u j ξ u ⊤ M(ξ, θ)u , C tot (ξ, θ) = sup u∈R p -{0} p i,j=1 u i [∂ 2 η(•, θ)/∂θ i ∂θ j ]u j ξ u ⊤ M(ξ, θ)u , ≤ C int (ξ, θ) + C par (ξ, θ) ,
with P θ the projector

(P θ f ) (x ′ ) = ∂η(x ′ , θ) ∂θ ⊤ M -1 (ξ, θ) X ∂η(x, θ) ∂θ f (x) ξ(dx) ,
and correspond to the original measures of nonlinearity of [START_REF] Bates | Relative curvature measures of nonlinearity[END_REF] for σ = 1, with an adaptation to the use of a design measure ξ instead of an exact design (x 1 , . . . , x N ). The connection with the curvature arrays of [START_REF] Bates | Relative curvature measures of nonlinearity[END_REF] is presented in (Pázman, 1993, Sect. 5.5); a procedure for their numerical computation is given in [START_REF] Bates | Relative curvature measures of nonlinearity[END_REF][START_REF] Ratkowsky | Nonlinear Regression Modelling[END_REF].

All computations are performed in Matlab on a bi-processor PC (2.5 Ghz) with 64 bits, equipped with 32 Go RAM. Classical optimal designs (D-, E-and c-optimality) are computed with the cutting-plane method, see Pronzato and Pázman (2013, Sect. 9.5.3); LP problems are solved with simplex algorithm; we use sequential quadratic programming for the local minimization of H(w k , θ) that yields θ (k+1) in ( 16)-(ii). Example 2. This example is artificial and constructed to illustrate the possible pitfall of using a local approach (here E-optimal design) for designing an experiment. The model response is given by

η(x, θ) = θ 1 {x} 1 + θ 3 1 (1 -{x} 1 ) + θ 2 {x} 2 + θ 2 2 (1 -{x} 2 ) , θ = (θ 1 , θ 2 ) ⊤ ,
with x ∈ X = [0, 1] 2 and {x} i denoting the i-th component of x. We consider local designs for θ 0 = (1/8, 1/8) ⊤ . One may notice that the set {∂η(x, θ)/∂θ θ 0 : x ∈ X } is the rectangle [3/64, 1] × [1/4, 1], so that optimal designs for any isotonic criterion function of the information matrix M(ξ) are supported on the vertices (0, 1), (1, 0) and (1, 1) of X . The classical Dand E-optimal designs are supported on three and two points respectively,

ξ * D,θ 0 ≃    0 1 1 0 1 1 0.4134 0.3184 0.2682    , ξ * E,θ 0 ≃    0 1 1 0 0.5113 0.4887    .
When only the design points x 1 = (0 1) ⊤ and x 2 = (1 0) ⊤ are used, the parameters are only locally estimable. Indeed, the equations in θ

′ η(x 1 , θ ′ ) = η(x 1 , θ) η(x 2 , θ ′ ) = η(x 2 , θ)
give not only the trivial solutions θ ′ 1 = θ 1 and θ ′ 2 = θ 2 but also θ ′ 1 and θ ′ 2 as roots of two univariate polynomials of the fifth degree (with coefficients depending on θ). Since these polynomials always admit at least one real root, at least one solution exists for θ ′ that is different from θ. In particular, the vector θ 0 ′ = (-0.9760 , 1.0567) ⊤ gives approximately the same values as θ 0 for the responses at x 1 and x 2 .

imsart-aos ver. 2013/03/06 file: Extended-optimality_AOS-REV.tex date: March 31, 2014

Direct calculations indicate that, for any θ, the maximum of η(•, θ)η(•, θ 0 ) 2 ξ with respect to ξ ∈ Ξ is reached for a measure supported on (0, 0), (0 , 1), (1, 0) and(1, 1). Also, the maximum of [η(x, θ) -η(x, θ 0 )] 2 with respect to x is attained on the same points. We can thus restrict our attention to the design space X = {(0, 0), (0, 1), (1, 0), (1, 1)}. We take Θ = [-3, 4] × [-2, 2] and use the algorithm of Sect. 3.4, with the grid G 0 of ( 16)-(iii) given by a random Latin hypercube design with 10 000 points in [0, 1] 2 renormalized to Θ (see, e.g., [START_REF] Tang | Orthogonal Array-Based Latin Hypercubes[END_REF]), to determine optimal designs for φ eE (•) and φ eG (•). When initialized with the uniform measure on the four points of X , and with ǫ = 10 -10 , the algorithm stops after 46 and 15 iterations, respectively requiring 0.67 s and 0.28 s in total, and gives the designs

ξ * eE,θ 0 ≃    0 0 0 1 1 1 0.32 0.197 0.483    , ξ * eG,θ 0 ≃    0 0 0 1 1 0 1 1 0.258 0.258 0.258 0.226    .
The performances of the designs ξ * D , ξ * E , ξ * eE and ξ * eG are given in Table 1. The values φ eE (ξ * E ) = φ eG (ξ * E ) = 0 indicate that E-optimal design is not suitable here, the model being only locally identifiable for ξ * E . The parametric, intrinsic and total measures of curvature at θ 0 (for σ 2 = 1) are also indicated in Table 1, see Pronzato and Pázman (2013, p. 223). Notice that the values of these curvature at θ 0 do not reveal any particular difficulty concerning ξ * E , but that the lack of identifiability for this design is pointed out by the extended optimality criteria.

This example is very particular and situations where the model is locally, but not globally, identifiable are much more common: in that case, (2) is only satisfied locally, for θ ′ in a neighborhood of θ, and one may refer, e.g., to Walter (1987); [START_REF] Walter | Identifiabilities and nonlinearities[END_REF] for a precise definition and examples. The lack of global identifiability would then not be detected by classical optimal design, but the maximum of φ eE (•) and φ eG (•) would be zero for Θ large enough, showing that the model is not globally identifiable.

Example 3. Consider the regression model (one-compartment with firstorder absorption input) used in [START_REF] Atkinson | Optimal experimental designs for properties of a compartmental model[END_REF], 3,6] and use the algorithm of Sect. 3.4 to compute an optimal design for φ eE (•); the grid G 0 of ( 16)-(iii) is taken as a random Latin hypercube design with 10 000 points in [0, 1] 3 renormalized to Θ. The number of iterations and computational time depend on ℓ, the number of elements of X . For instance, when X is the finite set {0.2, 0.4, 0.6, . . . , 24} with ℓ = 120, and the required precision ǫ equals 10 -10 , the algorithm initialized at the uniform measure on the three points 0.2, 1 and 23 converges after 42 iterations in about 26 s. By refining X iteratively around the support points of the current optimal design, after a few steps we obtain ξ * eE,θ 0 ≃ 0.1785 1.520 20.95 0.20 0.66 0.14 .

(19) η(x, θ) = θ 1 [exp(-θ 2 x) -exp(-θ 3 x)] , θ = (θ 1 , θ 2 , θ 3 ) ⊤ , x ∈ R + , imsart-
2; det 1/3 = φD(ξ) = {det[M(ξ, θ 0 )]} 1/3 , λmin = φE(ξ) = λmin[M(ξ, θ 0 )].
A similar approach is used below for the construction of optimal designs for φ ec (•) and in Example 4 for φ eE (•). The performances of the designs ξ * D , ξ * E and ξ * eE are indicated in Table 2. One may notice that the design ξ * eE is best or second best for φ D (•), φ E (•) and φ eE (•) among all locally optimal designs considered.

The intrinsic curvature is zero for ξ * D , ξ * E and ξ * eE (since they all have 3 = dim(θ) support points) and the parametric curvatures at θ 0 are rather small (the smallest one is for ξ * eE ). This explains that, the domain Θ being

(ξ) = λ min [M(ξ, θ 0 )].
Consider now the same three functions of interest as in [START_REF] Atkinson | Optimal experimental designs for properties of a compartmental model[END_REF]: g 1 (θ) is the area under the curve,

g 1 (θ) = ∞ 0 η(x, θ) dx = θ 1 (1/θ 2 -1/θ 3 ) ; g 2 (θ) is the time to maximum concentration, g 2 (θ) = log θ 3 -log θ 2 θ 3 -θ 2 ,
and g 3 (θ) is the maximum concentration

g 3 (θ) = η[g 2 (θ), θ] .
We shall write c i = c i (θ 0 ) = ∂g i (θ)/∂θ| θ 0 and denote ξ c i ,θ 0 the (locally) optimal design for g i (θ) which maximizes For each function g i , we restrict the search of a design ξ ec i optimal in the sense of the criterion φ ec (•) to design measures supported on the union of the supports of ξ * D,θ 0 , ξ * E,θ 0 and ξ * c i ,θ 0 . We then obtain the following designs .

φ c i (ξ; θ 0 ) = [c ⊤ i M -(ξ, θ 0 )c i ] -1 , for i = 1, 2, 3. The ξ * c i ,
ξ * ec 1 ,θ 0 ≃ 0.2327 1.389 23.36 9 • 10 -4 1.2 • 10 -2 0.9871 , ξ * ec 2 ,θ 0 ≃ 0.
The performances of ξ * c i and ξ * ec i , i = 1, . . . , 3, are indicated in Table 2, together with the curvature measures at θ 0 for ξ * ec i (which are nonsingular). For each function g i of interest, the design ξ * ec i performs slightly worse than ξ * c i in terms of c-optimality, but, contrarily to ξ * c i , it allows us to estimate the three parameters θ and guarantees good estimability properties for g i (θ) for all θ ∈ Θ. Notice that, apart from the c-optimality criteria φ c i (•), all criteria considered take the value 0 at the c-optimal designs ξ * c i . The construction of an optimal design for φ ec (•) thus forms an efficient method to circumvent the difficulties caused by singular c-optimal design in nonlinear models, see (Pronzato and Pázman, 2013, Chap. 3 & 5). One may also refer to [START_REF] Pronzato | On the regularization of singular c-optimal designs[END_REF] for alternative approaches for the regularization of singular c-optimal designs.

We conclude this example with a comparison with the average-optimal designs of [START_REF] Atkinson | Optimal experimental designs for properties of a compartmental model[END_REF], that aim at taking uncertainty on θ 0 into account. Consider a prior distribution π(•) on the two components of θ that intervene nonlinearly in η(x, θ), and let IE π {•} denote the expectation for π(•). [START_REF] Atkinson | Optimal experimental designs for properties of a compartmental model[END_REF] indicate that when π equals

π A uni- form on [θ 0 2 -0.01, θ 0 2 + 0.01] × [θ 0 3 -1, θ 0 3 + 1], the design that maximizes IE π {log det[M(ξ, θ)]} is ξ * AD-A ≃ 0.2288 1.4170 18.4513 1/3 1/3 1/3 ,
and the designs that minimize Their performances are indicated in the bottom part of Table 2. Notice that the average-optimal designs for the vague prior π B are supported on more than three points and thus allow model checking; this is the case too for the two designs ξ * ec 2 and ξ * ec 3 . However, contrary to average-optimal design, the number of support points of optimal designs for extended optimality criteria does not seem to increase with uncertainty measured by the size of Θ: for instance, when Θ = [θ 0 1 -5, θ 0 1 + 5] × [θ 0 2 -0.04, θ 0 2 + 0.04] × [θ 0 3 -4, θ 0 3 + 4], the optimal design for φ eE (•) is still supported on three points, approximately 0.1565, 1.552 and 19.73, receiving weights 0.268, 0.588 and 0.144 respectively.

IE π {c ⊤ i (θ)M -(ξ, θ)c i (θ)}, i = 1, 2, 3, are ξ * Ac 1 -A ≃ 0.
All average-optimal designs considered yield reasonably small curvatures at θ 0 , although larger than those for ξ * E,θ 0 and ξ * eE,θ 0 . The performances of ξ * AD-A and ξ * AD-B are close to those of ξ * D,θ 0 , and the most interesting features concern designs for estimation of functions of interest g i (θ). The designs ξ * c i ,θ 0 cannot be used if θ = θ 0 and are thus useless in practice. The average-optimal designs ξ * Ac i -B perform significantly worse than ξ * ec i ,θ 0 in terms of φ ec i (•) for i = 1, 2 and 3 and in terms of φ c i (•) for i = 2 and 3. On the other hand, the designs ξ * Ac i -A , constructed for the precise prior π A , perform significantly better than ξ * ec i ,θ 0 in terms of φ c i (•) for all i. Figure 3 presents φ c 3 (ξ; θ) as a function of θ, for the three designs ξ * Ac 3 -A (dashed line), ξ * Ac 3 -B (dash-dotted line) and ξ * ec 3 ,θ 0 (solid line), when θ 1 = θ 0 1 , θ 3 = θ 0 3 (left) and θ 1 = θ 0 1 , θ 2 = θ 0 2 (right). Note that the projection on the last two components of θ of the set Θ used for extended c-optimality is intermediate between the supports of π A and π B . Although average-optimal designs ξ * Ac i -A,B and extended-optimal designs ξ * ec i ,θ 0 pursue different objectives, the example indicates that they show some resemblance in terms of precision of estimation of g i (θ). The situation would be totally different in absence of global identifiability for g i (θ), a problem that would not be detected by average-optimal designs, see the discussion at the end of Example 2.

]; det 1/3 = φD(ξ) = {det[M(ξ, θ 0 )]} 1/3 , λmin = φE(ξ) = λmin[M(ξ, θ 0 )].
The optimal (maximum) values of the criteria are on the main diagonal and indicated in boldface. The bottom part of the table corresponds to the average-optimal designs of [START_REF] Atkinson | Optimal experimental designs for properties of a compartmental model[END_REF]. Example 4. For the same regression model ( 19), we change the value of θ 0 and the set Θ and take θ 0 = (0.773 , 0.214 , 2.09) ⊤ and Θ = [0, 5] × [0, 5] × [0, 5], the values used by [START_REF] Kieffer | Interval analysis for guaranteed nonlinear parameter estimation[END_REF]. With these values, from an investigation based on interval analysis, the authors report that for the 16-point design

ξ 0 = 1 2 • • • 16 1/16 1/16 • • • 1/16
and the observations y given in their Table 13.1, the LS criterion yη X (θ) 2 has a global minimizer (the value we have taken here for θ 0 ) and two other local minimizers in Θ. The D-and E-optimal designs for θ 0 are now given by ξ * D,θ 0 ≃ 0.42 1.82 6.80 1/3 1/3 1/3 , ξ * E,θ 0 ≃ 0.29 1.83 9.0 0.4424 0.3318 0.2258 .

Using the same approach as in Example 3, with the grid G 0 of ( 16)-(iii) obtained from a random Latin hypercube design with 10 000 points in Θ, we obtain ξ * eE,θ 0 ≃ 0.38 2.26 7.91 0.314 0.226 0.460 .

To compute an optimal design for φ eG (•), we consider the design space X = {0, 0.1, 0.2, . . . , 16} (with 161 points) and use the algorithm of Sect. 3.4 with the grid G 0 of ( 16)-(iii) taken as a random Latin hypercube design 7. Further extensions and developments.

3 ; det 1/3 = φD(ξ) = {det[M(ξ, θ 0 )]} 1/3 , λmin = φE(ξ) = λmin[M(ξ, θ 0 )].
7.1. An extra tuning parameter for a smooth transition to usual design criteria. The criterion φ eE (ξ; θ 0 ) can be written as (20)

φ eE (ξ; θ 0 ) = max{α ∈ R : η(•, θ) -η(•, θ 0 ) 2 ξ ≥ α θ -θ 0 2 , for all θ ∈ Θ} .
Instead of giving the same importance to all θ whatever their distance to θ 0 , one may wish to introduce a saturation and reduce the importance given to those θ very far from θ 0 , that is, consider for some K ≥ 0. Equivalently, φ eE|K (ξ; θ 0 ) = min θ∈Θ H E|K (ξ, θ), with

φ eE|K (ξ; θ 0 ) = max α ∈ R : η(•, θ) -η(•, θ 0 ) 2 ξ ≥ α θ -θ 0 2 1 + K θ -θ 0 2 ,
H E|K (ξ, θ) = η(•, θ) -η(•, θ 0 ) 2 ξ K + 1 θ -θ 0 2 .
As in Sect. 3.1, we obtain φ eE|K (ξ) = λ min [M(ξ)] in a linear model and, for a nonlinear model with Θ = B(θ 0 , ρ), lim ρ→0 φ eE|K (ξ; θ 0 ) = λ min [M(ξ, θ 0 )] for any K ≥ 0. Moreover, in a nonlinear model with no overlapping φ eE|K (ξ; θ 0 ) can be made arbitrarily close to λ min [M(ξ, θ 0 )] by choosing K large enough, whereas choosing K not too large ensures some protection against η X (θ)η X (θ 0 ) being small for some θ far from θ 0 . Also, properties of φ eE (•) such as concavity, positive homogeneity, existence of directional derivatives, see Sect. 3.2, remain valid for φ eE|K (•), for any K ≥ 0. The maximization of φ eE|K (•) forms a LP problem when both X and Θ are finite, see Sect. 3.3, and a relaxation procedure (cutting-plane method) can be used when Θ is a compact subset of R p , see Sect. 3.4.

A similar approach can be used with extended c-and G-optimality, which gives φ ec|K (ξ) = min θ ′ ∈Θ H c|K (ξ, θ ′ ) with

H c|K (ξ, θ) = η(•, θ) -η(•, θ 0 ) 2 ξ K + 1 |g(θ) -g(θ 0 )| 2 , and φ eG|K (ξ) = min θ∈Θ η(•, θ) -η(•, θ 0 ) 2 ξ K + 1 max x∈X [η(x, θ) -η(x, θ 0 )] 2 ,
for K a positive constant.

7.2. Worst-case extended optimality criteria. The criterion defined by φ M eE (ξ) = min θ 0 ∈Θ φ eE (ξ; θ 0 ) = min (θ,θ 0 )∈Θ×Θ H E (ξ, θ; θ 0 ) , see ( 6), (5), accounts for the global behavior of η(•, θ) for θ ∈ Θ and obliterates the dependence on θ 0 that is present in φ eE (ξ; θ 0 ). The situation is similar to that in Sect. 3, excepted that we consider now the minimum of H E with respect to two vectors θ and θ 0 in Θ × Θ. All the developments in Sect. 3 obviously remain valid (concavity, existence of directional derivative, etc.), including the algorithmic solutions of Sects. 3.3 and 3.4. The same is true for the worst-case versions of φ ec (•) and φ eG (•), respectively defined by φ M ec (ξ) = min (θ,θ 0 )∈Θ×Θ H c (ξ, θ; θ 0 ), see (17), and by φ M eG (ξ) = min (θ,θ 0 )∈Θ×Θ { η(•, θ) -η(•, θ 0 ) 2 ξ / max x∈X η(x, θ) -η(x, θ 0 ) 2 } , and for the worst-case versions of the extensions of previous section that include an additional tuning parameter K.

Note that the criterion φ M eE (•) may direct attention to a particularly pessimistic situation. Indeed, for Θ a compact set with nonempty interior and µ the Lebesgue measure on Θ, one may have min θ 0 ∈Θ φ eE (ξ; θ 0 ) = 0 for all designs ξ although µ{θ 0 ∈ Θ : φ eE (ξ ′ ; θ 0 ) > 0} = 1 for some design ξ ′ . This corresponds to a situation where the model is structurally identifiable, in the sense that the property (2) is generic but is possibly false for θ in a subset of zero measure, see, e.g., Walter (1987)). Example 2 gives an illustration.

Example 2 (continued). When the three polynomial equations θ ′ 1 -θ ′

1 3 = θ 1 -θ 3 1 , θ ′ 2 -θ ′ 2 2 = θ 2 -θ 2 2 , θ ′ 1 3 + θ ′ 2 2 = θ 3 1 + θ 2
2 are satisfied, then η(x, θ ′ ) = η(x, θ) for all x. Since these equations have solutions θ ′ = θ in Θ × Θ, φ M eE (ξ) = 0 for all ξ ∈ Ξ. On the other hand, max ξ∈Ξ φ eE (ξ; θ 0 ) > 0 w.p.1 when θ 0 is randomly drawn with a probability measure having a density with respect to the Lebesgue measure on Θ.

In a less pessimistic version of worst-case extended E-optimality we may thus consider a finite set Θ 0 ⊂ Θ for θ 0 , obtained for instance by random sampling in Θ, and maximize min θ 0 ∈Θ 0 φ eE (ξ; θ 0 ). 8. Conclusions. Two essential ideas have been presented. First, classical optimality criteria can be extended in a mathematically consistent way to criteria that preserve a nonlinear model against overlapping, and at the same time retain the main features of classical criteria, especially concavity. Moreover, they coincide with their classical counterpart for linear models. Second, designs that are nearly optimal for those extended criteria can be obtained by standard linear programming solvers, supposing that the approximation of the feasible parameter space Θ by a finite set is acceptable. A relaxation method, equivalent to the cutting-plane algorithm, can be used when Θ is a compact set with nonempty interior. Linear constraints on the design can easily be taken into account. As a by-product, this also provides simple algorithmic procedures for the determination of E-, c-or G-optimal designs in linear models with linear cost constraints.

As it is usually the case for optimal design in nonlinear models, the extended-optimality criteria are local and depend on a guessed value θ 0 for the model parameters. However, the construction of a globalized, worst-case, version enjoying the same properties is straightforward (Sect. 7.2).

Finally, we recommend the following general procedure for optimal design in nonlinear regression. (i) Choose a parameter space Θ correspondimsart-aos ver. 2013/03/06 file: Extended-optimality_AOS-REV.tex date: March 31, 2014 ing to the domain of interest for θ, select (e.g., randomly) a finite subset Θ 0 in the interior of Θ; (ii) for each θ 0 in Θ 0 compute an optimal design ξ * eE,θ 0 maximizing φ eE (ξ; θ 0 ) and a E-optimal design ξ * E,θ 0 maximizing φ E (ξ; θ 0 ) = λ min M(ξ, θ 0 ); (iii) if φ eE (ξ * eE,θ 0 ; θ 0 ) is close enough to φ E (ξ * E,θ 0 ; θ 0 ) for all θ 0 in Θ 0 , one may consider that the risk of overlapping, or lack of identifiability in Θ, is weak and classical optimal design that focuses on the precision of estimation can be used; otherwise, a design that maximizes min θ 0 ∈Θ 0 φ eE (ξ; θ 0 ) should be preferred. When the extended Goptimality criterion φ eG (•; θ 0 ) is substituted for φ eE (•; θ 0 ), the comparison in (iii) should be between φ eG (ξ * eG,θ 0 ; θ 0 ) and 1/dim(θ), see Sect. 5. Extended c-optimality can be used when one is interested in estimating a (nonlinear) function of θ, the comparison in (iii) should then be with c-optimality.

Fig 1 .

 1 Fig 1. Expectation surface Sη for θ ∈ Θ = [0, 1], r = 1 and u = umax = 7 π/4.

Fig 2 .

 2 Fig 2. HE(νu, θ) (left) and ηX (θ) -ηX (θ 0 ) (right) as functions of θ ∈ Θ = [0, 1] for r = 1, u = 2 ( dotted line), u = umax = 7π/4 ( dashed line) and u = u * = π ( solid line).

φ

  eE (ξ; θ 0 ) = λ min [M(ξ, θ 0 )] .imsart-aos ver. 2013/03/06 file: Extended-optimality_AOS-REV.tex date: March 31, 2014

  imsart-aos ver. 2013/03/06 file: Extended-optimality_AOS-REV.tex date: March 31, 2014

  imsart-aos ver. 2013/03/06 file: Extended-optimality_AOS-REV.tex date: March 31, 2014

4. 1 .

 1 Definition and properties. Consider the case where one wants to estimate a scalar function of θ, denoted by g(θ), possibly nonlinear. We assume that c = c(θ) = ∂g(θ) ∂θ θ=θ 0 = 0 . imsart-aos ver. 2013/03/06 file: Extended-optimality_AOS-REV.tex date: March 31

  imsart-aos ver. 2013/03/06 file: Extended-optimality_AOS-REV.tex date: March 31, 2014

  θ 0 are singular and are approximately given by

Fig 3 .

 3 Fig 3. φc 3 (ξ; θ) as a function of θ, for θ1 = θ 0 1 and θ3 = θ 0 3 ( left) and θ1 = θ 0 1 , θ2 = θ 0 2 ( right); ξ = ξ * Ac 3 -A in dashed line, ξ = ξ * Ac 3 -B in dash-dotted line and ξ = ξ * ec 3 ,θ 0 in solid line.

Table 1

 1 aos ver. 2013/03/06 file: Extended-optimality_AOS-REV.tex date: March 31, 2014 Performances of designs ξ * D , ξ * E , ξ * eE and ξ * eG and curvature measures at θ 0 in Example

  The optimal (maximum) values of the criteria are indicated in boldface.

	ξ	det 1/3 λmin	φeE	φeG	Cpar Cint	Ctot
	ξ * D ξ * E ξ * eE ξ * eG	0.652 0.273 0.625 0.367 0.453 8.45 • 10 -2 8.78 • 10 -3 9.74 • 10 -2 3.33 3.16 • 10 -3 0.108 1.10 0 0 1.19 0.540 0.195 5.68 • 10 -3 0.340 1.33	0.541 1.22 0 1.19 2.69 4.28 1.26 1.83
	with nominal parameters θ 0 = (21.80 , 0.05884 , 4.298) ⊤ . The D-and E-
	optimal designs for θ 0 are respectively given by
		ξ * D,θ 0 ≃	0.229 1.389 18.42 1/3 1/3 1/3	,
		ξ * E,θ 0 ≃	0.170 1.398 23.36 0.199 0.662 0.139	,

see

[START_REF] Atkinson | Optimal experimental designs for properties of a compartmental model[END_REF]

.

We take Θ as the rectangular region [16, 27]×[0.03, 0.08]×[

  imsart-aos ver. 2013/03/06 file: Extended-optimality_AOS-REV.tex date: March 31, 2014 not too large, the values of φ eE (ξ) do not differ very much from those of φ E

Table 2 .

 2 Performances of different designs and curvature measures at θ 0 for the model (19) with θ 0 = (21.80 , 0.05884 , 4.298) ⊤ and Θ = [16, 27] × [0.03, 0.08] × [3, 6

	imsart-aos ver. 2013/03/06 file: Extended-optimality_AOS-REV.tex date: March 31, 2014

Table 3

 3 Performances of different designs and curvature measures at θ 0 for the model (19) with θ 0 = (0.773 , 0.214 , 2.09) ⊤ and Θ =[0, 5]

  The optimal (maximum) values of the criteria are indicated in boldface.with 10 5 points. The same design space is used to evaluate φ eG (•) for the four designs above. For ǫ = 10 -10 , the algorithm initialized at the uniform measure on X converges after 34 iterations in about 52 s and gives The large intrinsic curvature for ξ 0 , associated with the small values of φ eE (ξ 0 ) and φ eG (ξ 0 ), explains the presence of local minimizers for the LS criterion and thus the possible difficulties for the estimation of θ. The values of φ eE (•) and φ eG (•) reported in the table indicate that ξ *

	ξ	det 1/3	λmin	φeE	φeG	Cpar	Cint	Ctot
	ξ0	1.85 • 10 -2	1.92 • 10 -4	2.28 • 10 -5	5.66 • 10 -3 180.7 15.73 181.3
	ξ * D ξ * E ξ * eE ξ * eG	5.19 • 10 -2 1.69 • 10 -3 4.51 • 10 -2 2.04 • 10 -3 1.32 • 10 -4 2.64 • 10 -4 4.73 • 10 -2 1.53 • 10 -3 2.92 • 10 -4 0.114 6.70 • 10 -2 58.0 7.95 • 10 -2 50.7 54.6 4.11 • 10 -2 1.31 • 10 -3 1.69 • 10 -4 0.244 69.7	0 0 0 10.7	58.0 50.7 54.6 69.9
			ξ * eG,θ 0 ≃	0.4 0.278 0.258 0.244 0.22 1.9 5.3 16	.
	ξ * eE or ξ * eG would have caused less difficulties.			D , ξ * E ,

The performances and curvature measures at θ 0 of ξ 0 , ξ * D , ξ * E , ξ * eE and ξ * eG are given in Table

3

.

imsart-aos ver. 2013/03/06 file: Extended-optimality_AOS-REV.tex date: March 31, 2014

Acknowledgements. The authors thank the referees for useful comments that helped to significantly improve the paper.

The work of Andrej Pázman was supported by the VEGA grant No.