Hakim Sidahmed

Elena Prokofyeva

Matthew B Blaschko

Discovering Predictors of Mental Health Service Utilization with k-support Regularized Logistic Regression

Keywords: variable selection, discarding variables, k-support regularized logistic regression, epidemiological data

Many epidemiological studies are undertaken with a use of large epidemiological databases, which involves the simultaneous evaluation of a large number of variables. Epidemiologists face a number of problems when dealing with large data sets: multicolinearity (when variables are correlated to each other), confounding factors (when risk factor is correlated with both exposure and outcome variable), and interactions (when the direction or magnitude of an association between two variables differs due to the effect of a third variable). Correct variable selection helps to address these issues and helps to obtain unbiased results. Selection of relevant variables is a complicated and a time consuming task. Flawed variable selection methods still prevail in the scientific literature; there is a need to demonstrate the usability of new algorithms using real data. In this paper we propose to use a novel machine learning method, k-support regularized logistic regression, for discovering predictors of mental health service utilization in the National Epidemiologic Survey for Alcohol and Related Conditions (NESARC). We show that k-support regularized logistic regression yields better prediction accuracy than 1 or 2 regularized logistic regression as well as several baseline methods on this task, and we qualitatively evaluate the top weighted variates. The selected variables are supported by related epidemiological research, and give important cues for public policy.

Introduction

Statistical analysis, such as logistic regression, is regularly applied to problems in epidemiology, and other problems where associations must be tested in the context of confounding factors. Often, such models are applied to data consisting of few variables hand-selected from a large, unstructured database. The validity of such studies may therefore depend on scientists carefully identifying all possible confounding factors. Furthermore, important related factors to a condition of interest may be excluded from analysis erroneously. It is therefore of interest to consider strategies that may include a large amount of available data, and select relevant variables automatically. This paper proposes a novel method for this kind of analysis, k-support regularized logistic regression, and we demonstrate its applicability on the challenging task of identifying predictors of mental health service utilization for people with substance abuse.

One problem that arises when using statistical learning on large databases is that spurious correlation with the signal is more likely with a large number of variables in a given analysis. A strategy to counteract this problem is to make use of a sparsity regularizer, such as the lasso, or 1 penalty [START_REF] Tibshirani | Regression shrinkage and selection via the lasso[END_REF][START_REF] Sha | Multiplicative updates for l1-regularized linear and logistic regression[END_REF][START_REF] Yuan | An improved GLMNET for l1-regularized logistic regression[END_REF]. While the lasso is beneficial in setting many coefficients to zero, and therefore reducing the risk that the model rely on irrelevant variables, a problem arises when there are correlated signals [START_REF] Manninen | Leukemia prediction using sparse logistic regression[END_REF]. If two or more variables are discriminative with respect to the target of interest, the lasso has a tendency to select only one. This can be detremental in epidemiological studies for several reasons: (i) if the selected variable is a confounding factor, the causal factor will be suppressed, and (ii) if the selected variable is a causal factor, confounding factors may be suppressed and the effect of the causal factor will be over-estimated. The solution to this problem is to use a correlated sparsity regularizer.

Correlated sparsity is a property achieved by a subset of structured sparsity regularizers [START_REF] Zou | Regularization and variable selection via the elastic net[END_REF][START_REF] Huang | Learning with structured sparsity[END_REF][START_REF] Jenatton | Structured variable selection with sparsity-inducing norms[END_REF]. This family of regularizers incorporates those that favor sparse coefficient vectors, but those that have specific patterns of sparsity. In this work, we consider specifically the k-support norm [START_REF] Argyriou | Sparse prediction with the k-support norm[END_REF], which achieves a lower degree of sparsity than the lasso, but favors configurations of the model in which variables with correlated signals are allowed to be activated simulatneously. We consider the novel application of k-support norm regularization to logistic regression. The ksupport norm is related mathematically to the elastic net [START_REF] Zou | Regularization and variable selection via the elastic net[END_REF], but has a more principled derivation as the convexification of a norm that applies an 2 penalty to a subset of k-selected variables.

Current Practices for Variable Selection in Epidemiological Research

The main aim of epidemiological studies is to identify and to quantify risk factors of diseases. Modern epidemiological research is often based on a large data sets such as the National Epidemiologic Survey for Alcohol and Related Conditions (NESARC) [START_REF] Chen | Alcohol use and alcohol use disorders in the United States: Main findings from the 2001-2002 National Epidemiologic Survey on Alcohol and Related Conditions (NESARC)[END_REF][START_REF] Prokofyeva | The role of family history in mental health service utilization for major depression[END_REF][START_REF] Melchior | Treatment for illegal drug use disorders: the role of comorbid mood and anxiety disorders[END_REF], the National Health and Nutritional Examination Survey (NHANES) [START_REF] Chakravorty | Daytime sleepiness: Associations with alcohol use and sleep duration in Americans[END_REF]; the CONSTANCES cohort [START_REF] Zins | The CONSTANCES cohort: an open epidemiological laboratory[END_REF], and the Whitehall cohort study [START_REF] Marmot | Health inequalities among British civil servants: the Whitehall II study[END_REF][START_REF] Hinnouho | Metabolically healthy obesity and the risk of cardiovascular disease and type 2 diabetes: The Whitehall II cohort study[END_REF], etc. Problems that scientists face when dealing with such data sets include: multicolinearity (when variables are correlated to each other), confounding factors (when a risk factor is correlated with both exposure and outcome variable), interactions (when the direction or magnitude of an association between two variables differs due to the effect of a third variable), the sample size (the study can be too large with meaningful associations being declared or too small to detect important associations), and the number of factors being studied (the higher the number of factors the higher is the probability to find interactions due to chance alone) [START_REF] Dohoo | An overview of techniques for dealing with large numbers of independent variables in epidemiologic studies[END_REF]. Correct variable selection helps to overcome these problems and to obtain unbiased results. In particular, it is used for confounder control in etiologic research [START_REF] Walter | Variable selection: Current practice in epidemiological studies[END_REF] and for unbiased estimation of probabilities in prediction research [START_REF] Greenland | Invited commentary: Variable selection versus shrinkage in the control of multiple confounders[END_REF].

It is often difficult to determine which variables are important and should be included in the analysis and which should be disregarded. Selection of relevant variables depends on number of factors: the research question of the study, study design, and the sample size [START_REF] Walter | Variable selection: Current practice in epidemiological studies[END_REF]. Prior scientific knowledge was earlier used as the main criteria to identify covariates that should be included in the analysis, but it is not always available for all research questions asked [START_REF] Greenland | Invited commentary: Variable selection versus shrinkage in the control of multiple confounders[END_REF][START_REF] Hernán | Causal knowledge as a prerequisite for confounding evaluation: An application to birth defects epidemiology[END_REF][START_REF] Walter | Variable selection: Current practice in epidemiological studies[END_REF]. The first step of any epidemiological study planning is the review of the prior scientific knowledge, in particular, other epidemiological studies that were performed on a similar subject, but also those for related pathologies. This is done to identify different factors that can influence the outcome of the study along with the risk factor of interest. Although this strategy of identifying covariates is an obligatory step for planning for any epidemiological study, it may not be suitable for all studies, for example, for those with a research question that has not previously been studied.

A number of statistical techniques were developed to enable the selection of covariates based on relations of the data under study: change in the effect estimate, stepwise selection, shrinkage and penalized regression, and other techniques. Recent research on contemporary epidemiological analysis showed that effect estimate change and stepwise selection techniques represent the majority of the all methods used for variable selection [START_REF] Walter | Variable selection: Current practice in epidemiological studies[END_REF]. Despite the fact that they allow the automatic selection of important predictor variables for inclusion in the model through the forward, backward or stepwise selection process, they have number of limitations. Firstly, these methods assume independence of all variables studied; secondly, there are many different ways to decide which model suits this data the best, and thirdly these methods are time consuming and often require greater input from the researcher [START_REF] Dohoo | An overview of techniques for dealing with large numbers of independent variables in epidemiologic studies[END_REF]. Indeed, for experiments on the scale considered here, forward and backward selection is computationally infeasible (see Section 3.3.1). In the light of the dominance of variable selection methods that are criticized by leading epidemiologists as flawed, there is a need to demonstrate the usability of new algorithms in real data instead of simulation studies [START_REF] Walter | Variable selection: Current practice in epidemiological studies[END_REF].

Alcohol abuse and dependence remain one of the main public health problems. It is closely linked to car crashes [START_REF] Chou | The prevalence of drinking and driving in the United States, 2001-2002: Results from the national epidemiological survey on alcohol and related conditions[END_REF], domestic violence [START_REF] Caetano | Intimate partner violence, dependence symptoms and social consequences from drinking among white, black and hispanic couples in the United States[END_REF], fetal alcohol syndrome [START_REF] Lemoine | Children of alcoholic parents-observed anomalies: Discussion of 127 cases[END_REF], neurophysiological impairment [START_REF] Bates | Neurocognitive impairment associated with alcohol use disorders: Implications for treatment[END_REF], and psychiatric comorbidities [START_REF] Grant | Prevalence and co-occurrence of substance use disorders and independentmood and anxiety disorders: Results from the national epidemiologic survey on alcohol and related conditions[END_REF]. Despite the obvious treatment benefits such as reduction of the severity and harm associated with alcohol use disorder, only 7% of individuals with lifetime alcohol use disorder receive treatment [START_REF] Hasin | Prevalence, correlates, disability, and comorbidity of DSM-IV alcohol abuse and dependence in the United States: Results from the national epidemiologic survey on alcohol and related conditions[END_REF]. Therefore it is essential to study predictors of umet need for treatment for alcohol use disorder.

The aim of this article is to demonstrate the application of k-support regularized logistic regression to discover predictors of treatment or the perceived unmet need for treatment among respondents did not receive treatment for alcohol use disorder.

k-support Regularized Logistic Regression

In the sequel we will use the notation β ∈ R d to specify a coefficient vector of a linear discriminant function, λ ≥ 0 to specify a scalar regularization parameter, Ω : R d → R + a regularizer for a linear function class. We assume a training set of data {(x i , y i )} 1≤i≤n ∈ (X × Y) n which are encoded in a matrix X ∈ R d×n and a matrix Y ∈ {-1, +1} c×n where c is the number of classes. We assume i.i.d. sampling of this training set. The ith column of X is equal to x i and the ith column of Y encodes the class of y i by a binary vector with exactly one positive entry in the place of the corresponding class. For the special case of binary classification, i.e. c = 2, we will use the vector y ∈ {-1, +1} n to encode the class labels.

Regularized Logistic Regression

In the case of binary classification, logistic regression can be written as a regularized risk minimization with logistic loss. In this case, the objective to be minimized is

β = arg min β λΩ(β) + f log (β, X, y) (1)
with

f log (β, X, y) = n i=1 log 1 + e -y i β,x i (2)
where λ is a scalar parameter controlling the degree of regularization by the regularizer Ω and •, • is the canonical inner product in R d . The minimization of logistic loss is equivalent to classicial logistic regression with two classes in the case that it is unregularized [START_REF] Hastie | The Elements of Statistical Learning: Data Mining, Inference, and Prediction[END_REF].

In the multiclass setting, we specify one coefficient vector per class. With |Y| = c classes, we write our joint probability

p(β)p(Y|X; β) = 1 Z λΩ e -λΩ(β) • n i=1 e β y i ,x i c l=1 e β l ,x i (3)
where Z λΩ is a constant normalizing e -λΩ(β) as a distribution. Taking negative logarithms we arrive at

-log p(β)p(Y|X; β) = log Z λΩ + λΩ(β) + n i=1        -β y i , x i + log c l=1 e β l ,x i        (4)
We ignore the term log Z λΩ as there is no dependence on β, and we arrive at the loss function

f log (β, X, Y) = n i=1        -β y i , x i + log c l=1 e β l ,x i        (5)
where we have abused the notation f log . When the matrix Y is used as an argument to f log , we signify that we are using the definition in Equation ( 5), while the use of the vector y indicates that we are in the binary setting and using Equation (2)

In the multi-class setting here we have written our joint probability slightly differently than the classical presentation of logistic regression [START_REF] Hastie | The Elements of Statistical Learning: Data Mining, Inference, and Prediction[END_REF]Equation (4.18)]. In that presentation, one of the classes had its class energy arbitrarily constrained to 1 as the problem is overparametrized when performing maximum likelihood learning on the probability simplex:

p(y = i|x; β) = e β i ,x 1 + c-1 l=1 e β l ,x ∀i c (6) p(y = c|x; β) = 1 1 + c-1 l=1 e β l ,x (7)
In the case of maximum a posteriori (MAP) estimation, we would like our structured-sparsity prior over functions to be symmetric with respect to the classes so that each will be regularized equally:

p(y = i|x; β) = e β i ,x c l=1 e β l ,x ∀i (8)
This formulation is ill posed in the maximum likelihood setting, but the use of a regularizer yields a well-posed symmetric problem for MAP estimation.

k-support Regularization

We propose the use of a correlated sparsity regularizer, the k-support norm [START_REF] Argyriou | Sparse prediction with the k-support norm[END_REF]. The k-support norm is the gauge function associated with the convex set

conv{β | β 0 ≤ k, β 2 ≤ 1}. ( 9
)
This set is the convex hull of all d k axis aligned unit 2 balls of dimensionality k. k is a parameter that will be chosen by model selection or prior information.

The k-support norm can be computed as [START_REF] Prokofyeva | The role of family history in mental health service utilization for major depression[END_REF] where |β| ↓ i is the ith largest element of the vector and r is the unique integer in {0, . . . , k -1} satisfying

β sp k =          k-r-1 i=1 (|β| ↓ i ) 2 + 1 r + 1         d i=k-r |β| ↓ i         2          1 2
|β| ↓ k-r-1 > 1 r + 1 d i=k-r |β| ↓ i ≥ |β| ↓ k-r . (11)
One may observe that the k-support norm computes a weighted sum of the 2 norm on the largest components of the vector (leftmost sum in Equation ( 10)) and the 1 norm on the smallest components of the vector (rightmost sum in Equation ( 10)). The integer r is dependent on k, but can increase the number of non-zero coefficients to be larger than k depending on where on the k-support ball the vector β is located.

Relationship to Other Methods

Relationship to 1 and 2 Regularization

In the case that k = 1 the k-support norm is exactly equivalent to the 1 norm. In the case that k = d, where β ∈ R d , the k-support norm is equivalent to the 2 norm.

We note that for an objective min

β λ β sp k + f (β, X, y) (12)
with some loss function

f (•, •, •), when k = d, this is equivalent to min β λ β 2 + f (β, X, y) (13)
rather than the familiar squared 2 regularizer. However, for any λ there exists some λ such that arg min

β λ β 2 + f (β, X, y) = arg min β λ β 2 2 + f (β, X, y). (14)
This can be easily seen by noting that the objectives are the Lagrangians of constrained minimization problems that minimize f subject to the equivalent constraints β 2 ≤ B and β 2 2 ≤ B 2 , respectively, for some B ∈ R + . We note that 1 and 2 reglarized logistic regression are precisely those methods considered in [START_REF] Ng | Feature selection, l 1 vs. l 2 regularization, and rotational invariance[END_REF] and that our proposed method therefore specializes to these algorithms when k = 1 or k = d, respectively.

Relationship to Group Lasso

A prominent family of structured sparsity regularizers are those based on the group lasso with overlaps [START_REF] Jacob | Group lasso with overlap and graph lasso[END_REF][START_REF] Jenatton | Structured variable selection with sparsity-inducing norms[END_REF]. These methods generalize the 1 regularizer to an 1 norm of a vector of 2 norms computed on potentially overlapping groups of variables. As such, one may view the k-support norm regularization as a special case of the group lasso in which groups are defined over all overlapping d k groups of size k. The group lasso with non-overlapping groups for logistic regression has been proposed in [START_REF] Meier | The group lasso for logistic regression[END_REF]. The method proposed here can be viewed as a special case of group lasso with overlapping groups and logistic loss. The main benefits of this special case are that the k-support norm has a beneficial interpretation in terms of correlated sparsity [START_REF] Argyriou | Sparse prediction with the k-support norm[END_REF] and can be computed in O(d log d) time, while optimization strategies for arbitrarily defined groups will scale with the exponential number of groups.

Optimization

Our method for optimizing k-support regularized logistic regression essentiallly follows the strategy suggested for k-support regularized least squares in [START_REF] Argyriou | Sparse prediction with the k-support norm[END_REF]. The strategy is based on Nesterov's accelerated method for the linear combination of smooth and non-smooth functions, where f log is the smooth component and • sp k non-smooth. Accelerated methods were first proposed for smooth functions [START_REF] Nesterov | A method of solving a convex programming problem with convergence rate O(1/k 2 )[END_REF][START_REF] Nesterov | On an approach to the construction of optimal methods of minimization of smooth convex functions[END_REF] and later applied to non-smooth [START_REF] Nesterov | Smooth minimization of non-smooth functions[END_REF] and composite functions [START_REF] Nesterov | Gradient methods for minimizing composite objective function[END_REF][START_REF] Beck | A fast iterative shrinkage-thresholding algorithm for linear inverse problems[END_REF]. An overview with a unified analysis is presented in [START_REF] Tseng | On accelerated proximal gradient methods for convex-concave optimization[END_REF].

The optimization strategy requires as inputs a function that gives oracle access to the smooth component, a function that gives the gradient of the smooth function, and a Lipschitz constant of the gradient:

f log (β, X, y) = n i=1 log 1 + e -y i β,x i (15) ∂ f log ∂β = - n i=1 e -y i β,x i 1 + e -y i β,x i y i x i (16)
L log = γ 4 ( 17
)
where the Lipschitz constant has a factor 1 4 from the Lipschitz constant of the sigmoid in

∂ f log
∂β , and γ is the largest eigenvalue of XX T . To minimize J(β) = λ β sp k + f log (β, X, y), convergence of accelerated methods is such that at the mth iteration

J β (m) -J (β * ) ≤ O 1 t min (m + 1) 2 (18)
where we have suppressed dimensionality dependent terms in the notation O(1), β (m) is the estimate of the coefficient vector after m steps of the optimization algorithm, β * is the optimal coefficient vector

β * := arg min β J(β) (19)
and

t min = min 1, C L log ( 20
)
for a method dependent constant 0 < C ≤ 1. The convergence therefore scales proportionately with the Lipschitz constant of the gradient of the loss L log . The Lipschitz constant of the gradient of the loss is largely dependent on the largest eigenvalue of XX T , and so the convergence matches the usual dependency on the condition number of the system matrix for second order optimization methods [START_REF] Boyd | Convex Optimization[END_REF]. The optimization for the multi-class variant follows

analogously. An open source implementation of k-support regularized logistic regression is available for download from https://github.com/blaschko/ksupport.

Results

The Epidemiological Data Set

The NESARC is a nationally representative study of the US population that was conducted among non-institutionalized adults (≥ 18 years of age) residing in households and group quarters designed and conducted by the National Institute We have applied k-support regularized logistic regression to the problem of predicting whether a given NESARC participant had an unmet need for alcohol addiction treatment. This condition was identified based on the subject's response to the question S2CQ4A in the NESARC survey. In the sequel, we will refer to this response as variable A. This convention (and that for other NESARC variables referred to in the text) is specified in Table 1.

Data Processing

It is primarily of interest to determine potential causes for not receiving treatment for alcoholism for those subjects that have ever had any alcohol dependence. This information is encoded in two variables in the NESARC database, one which encodes whether the respondent had any alcohol abuse or dependence in the last 12 months, or whether they had any alcohol or dependence prior to the last 12 months. We therefore filter respondents to only those that answered yes to variable B or variable C (cf. Table 1). Out of the original 43093 respondents, 4068 had suffered from alcoholism at some point prior to the administration of the survey.

Wave 1 of NESARC contains 2991 variables in total. By using a standard epidemiological approach based on prior scientific knowledge of the factors that may have an influence on help seeking for alcohol use disorder we were able to reduce our data set to 112 potentially relevant variables. This process consisted of performing a literature review of relevant variables included in previous studies and then including all variables in the relevant NESARC section, giving a substantially larger number of variables to the statistical learning process than previous studies. Finally, we removed variables corresponding to questions only relevant to people who do not receive treatment for alcohol dependence or abuse. Inclusion of these variables in the regression would effectively give the learning algorithm access to the label for a large number of respondents and render the output invalid. This pre-selection procedure is therefore essential to maintain the validity of the result, and has the added effect of making use of epidemiological prior knowledge. We then whiten the resulting variables to have zero mean and unit variance. In this way, the variable selection procedure encoded in the k-support regularization will not be biased towards any particular variables.

As is the case with surveys with a very large sample size and many questions, the data are incomplete. Many questions give the respondent the option to answer 'Unknown,' meaning that the respondent doesn't know the answer to the question (e.g. doesn't remember) and/or 'Blank,' when the respondent did not answer. These two cases were treated as distinct by encoding a binary variable indicating whether either event occurred.

Finally, we encode the categorical variables as follows: If a feature is binary (Yes(1)-No(0) answer), it is replaced by -1 (No) or +1 (Yes). If a feature is multivariate, it is expanded into several columns (corresponding to its number of classes) of 0's and 1's. If a feature takes values a, b, and c for instance, then it is expanded into 3 columns: a is replaced by vector [1, 0, 0] T , b by [0, 1, 0] T and c by [0, 0, 1] T . This encoding scheme results in a statistical learning problem with 314 dimensions.

The source code for data processing is available from https://github.com/hakimsd9/predUsingksup/.

Quantitative Experiments

We performed multiple experiments with random splits of the data into training (2,000 respondents), validation (1,000 respondents) and test (1,068 respondents) sets. For each experiment, there are two intrinsic parameters to choose: the regularization parameter λ and the value of k in the k-support norm. These parameters are selected using model selection, and the best performing configuration (as measured by area under the curve-AUC) on the validation set is used to train a discriminant function that is evaluated only once on the test set. We have used AUC as a selection criterion due to its better ability to distinguish between two methods when the accuracy is tied. Furthermore, recent theoretical and empirical results indicate that AUC optimization leads to good classification performance with tight regret bounds, further validating this choice of selection criterion [START_REF] Narasimhan | On the relationship between binary classification, bipartite ranking, and binary class probability estimation[END_REF].

Each experiment we run to apply logistic regression to the data set has to go through several iterations in order to choose the parameter k and the regularization parameter λ giving the best area under curve (AUC). We use a grid search over these variables to select the best configuration. The parameter k is allowed to take values in the set {2, 4, 8, 16, 32, 64, 128, 256}, while λ is in the set {10 -14 , . . . , 10 14 } in powers of 10. Both parameter ranges are therefore sampled logarithmically.

Comparison to baseline methods

We have first compared k-support regularized logistic regression to two baselines: ridge regression, and the selection of the single best variable. The ridge regression regularization parameter was selected from the same range of values of λ as for the logistic regression experiments. The single best variable was selected by finding the variable (or its negation) that gave the best AUC on the validation set. ROC curves for k-support regularized logistic regression, ridge regression, and the single best variable are shown in Figure 1.

We have performed statistical significance testing to formally validate the difference in performance between methods. First, performed 50 different trials to estimate the accuracy for each method. A Pearson χ 2 goodness of fit test was applied to evaluate whether the performances across trials were Gaussian distributed. This test showed that the data were not normally distributed (p ≈ 10 -3 ). Consequently, we have applied non-parametric significance tests here and in subsequent sections. We performed pairwise Wilcoxon signed rank tests, which showed significance for all pairs of variables (k-support vs. ridge regression p ≈ 2 × 10 -4 , k-support vs. single best variable p ≈ 5 × 10 -7 ). In addition to pairwise significance testing, we have also applied a Kruskal-Wallis test on all three methods simultaneously, showing a significant difference with p ≈ 3 × 10 -15 .. Forward and backward selection methods, although interesting potential baseline methods, are not feasible to apply in the context of this problem. This is because a model must be trained d(d+1) 2 -1 times, where d is the number of dimensions to be considered. As our problem contains 314 dimensions (see Section 3.2), this would require training a model 49,454 times, which is not possible in a reasonable amount of time. By contrast, sparsity regularization methods need only be trained once, resulting in feasible computation even for problems with a large number of variables. Recent theoretical results indicate that the conditions for forward selection to provide good variable selection are very close to those for 1 regularization [START_REF] Zhang | On the consistency of feature selection using greedy least squares regression[END_REF]. Consequently, we provide a comparison to 1 regularized logistic regression in the next section, which functions as a computationally tractable surrogate for forward or backward selection methods.

Comparison to 1 and 2 regularization

In order to compare the quality of the predictions obtained with the k-support regularizer to 1 and 2 regularizers, we perform Wilcoxon signed-rank tests as well as a Kruskal-Wallis test as in Section 3.3.1. The pairwise results are presented in Figure 2. The Kruskal-Wallis test showed significance with p ≈ 2 × 10 -4 . To illustrate the performance of the k-support regularized classifier, we plot its receiver operating characteristic curve (ROC curve). ROC curves for k-support regularization, 1 regularization, and 2 regularization are shown in Figure 3.

The variables with highest weight in the model learned by k-support regularization are summarized in Table 2. The obtained list of variables contained covariates that fell into five major groups listed in the order of relevance: 1) lifetime drug use disorder (covariates describing presents or absence of abuse/dependence on opioids, sedatives, tranquilizers, amphetamine, cocaine, inhalants/solvents, hallucinogens, cannabis, and heroin), 2) lifetime mental health problems (anxiety disorders: panic disorder with or without agoraphobia, agoraphobia with no history of panic disorder, generalized anxiety; mood disorders: manic disorder, dysthymia, hypomanic disorder, social phobia), and lifetime

              k-sup 1 2 k-sup 1 0.2059 1.9209e -06 1 0.2059 1 1.3536e -04 2 1.9209e -06 1.3536e -04 1                         k-sup 1 2
mean 0.7315 0.7265 0.6922 std 0.0255 0.0328 0.0098 3: ROC curves for k-support norm regularization, 1 regularization, and 2 regularization of logistic regression. k-support and 1 regularization perform substantially better than 2 regularization. k-support regularization performs slightly better than 1 regularization on average, but without statistical significance (see Figure 2). However, k-support regularization is to be favored due to its better handling of correlated variables (see Section 3.3.4). alcohol disorder, 3) family history of alcohol use disorder (full siblings, adoptive parents, natural mother, and natural grandmother) and family history of drug use disorders (natural mother, full siblings, and grandparent) and family history of mental health problems (natural mother, full sister), 4) type of treatment received for alcohol related problems (outpatient, inpatient, detoxification, rehabilitation, social services), 5) treatment received for general anxiety disorder

         

Convergence timing

We evaluate here the empirical performance of the optimization scheme described in Section 2.4 applied to ksupport regularized logistic regression. After fixing λ according to the best performing setting found in the previous, we have plotted the difference in primal objective over time for varying k (Figure 4). The special case of k = 1 corresponds to FISTA optimization [START_REF] Beck | A fast iterative shrinkage-thresholding algorithm for linear inverse problems[END_REF] of 1 regularized logistic regression, giving a strong state-of-the-art baseline method for the relative comparison of empirical convergence. We see that the convergence slows as k is increased, as has previously been noted by [START_REF] Argyriou | Sparse prediction with the k-support norm[END_REF], but does not exceed one order of magnitude over the 1 optimization even for high values of k. In practice, we are likely to be interested in low values of k, which exclude irrelevant variables and have given the best generalization performance in these experiments. All computation is feasible to perform in a short amount of time on a single core, and the convergence matches the theoretical bounds given in Equation (18).

Empirical correlation of selected variables

In order to test whether k-support regularized logistic regression is able to automatically select multiple correlated variables, we have computed the empirical correlation between selected variables and display this matrix in graphical form in Figure 5. The variable names correspond to the NESARC Codebook available from http: //pubs.niaaa.nih.gov/publications/AA70/AA70.htm. We first note that most pairwise correlations are close to zero, as is expected from a sparsity regularizer. However, for key subsets of related variables, there is high offdiagonal correlation. This gives strong support that the proposed correlated sparsity regularizer has the desired be- We use the same k-support code to compute 2 regularization (k = d), although in practice one could make use of optimization methods for twice-differentiable objective functions in this special case. There is a dependence between k and the speed of optimization, due to the longer computation required per iteration for larger values of k [START_REF] Argyriou | Sparse prediction with the k-support norm[END_REF]. All computation times are feasible even on a single core for standard commercial hardware, with convergence ranging from 1-35 seconds depending on the value of k. Figure 5: Empirical correlations between the top variables selected by k-support regularized logistic regression. Variable names are from the NESARC Codebook. We note that although there is a high degree of independence between variables, there are also significant correlations between some selected variables, a key benefit of the proposed method. The empirical behavior of the selected variables therefore matches the expectation of what should be achieved by a properly functioning correlated sparsity regularizer. See Section 3.3.4 for a discussion of these selected variables.

havior and is able to select multiple related variables. Furthermore, the variables for which there is a high degree of correlation indicate that the method has been able to automatically account for confounding variables. We observe high correlations between such pairs of selected variables as for example: dysthymia prior to the last 12 months (illness-induced and substance-induced ruled out) (DDYSSNISP12) and emergency room treatment because of dysthymia (S4CQ15B); hospitalization due to dysthymia (S4CQ15A) and emergency room visit due to dysthymia (S4CQ15B). The level of comorbidity of alcohol use disorders and mood disorders such as dysthymia is high: persons mood disorders have a 2-to 3-fold increased risk of alcohol use disorders [START_REF] De Graaf | Temporal sequencing of lifetime mood disorders in relation to comorbid anxiety and substance use disorders[END_REF][START_REF] Hasin | Prevalence, correlates, disability, and comorbidity of DSM-IV alcohol abuse and dependence in the United States: Results from the national epidemiologic survey on alcohol and related conditions[END_REF]. Secondary alcohol use disorders may result from self-medication of symptoms of dysthymia with alcohol [START_REF] Bolton | Self-medication of mood disorders with alcohol and drugs in the National Epidemiologic Survey on Alcohol and Related Conditions[END_REF], while dysthymia may be consequences of alcohol intoxication and/or withdrawal in primary alcohol use disorders [START_REF] Falk | Age of onset and temporal sequencing of lifetime DSM-IV alcohol use disorders relative to comorbid mood and anxiety disorders[END_REF]. Individuals who simultaneously suffer from alcohol use disorder and dysthymia may prefer to seek treatment for their dysthymic disorder, because of lesser stigma, easier access, and simpler reimbursement scheme. This will result in an increase of unmet need for treatment of their alcohol related problems [START_REF] Melchior | Treatment for illegal drug use disorders: the role of comorbid mood and anxiety disorders[END_REF]. Therefore, dysthymia can be considered as a known confounder.

Treatment in the emergency room because of dysthymia (S4CQ15B) may be considered as a potential cause of unmet need for alcohol treatment, due to the absence of an integrated treatment of comorbid psychiatric and addictive conditions [START_REF] Melchior | Treatment for illegal drug use disorders: the role of comorbid mood and anxiety disorders[END_REF]. Alcohol use disorder is often not treated by mental health services not specialized in addiction, therefore individuals with such comorbid conditions may have the impression that their alcohol use problem is not serious enough to require treatment or that adequate treatment is not available.

A second association is represented by two different types of treatment for dysthymia: hospitalization due to dysthymia (S4CQ15A) and emergency room visit due to dysthymia (S4CQ15B), in particular its more acute forms. Indeed, inclusion of these two variables instead of only one can help us to understand how types of treatment of dysthymia influence the unmet treatment for the alcohol use disorder. For the purpose of epidemiological analysis, these two variables will often be grouped together to represent a variable for treatment of more acute forms of dysthymia.

These observations confirm one of the key hypotheses of the application of k-support norm regularization in an epidemiological task: that confounders and other correlated variables can be automatically selected, strengthening the epidemiological analysis.

Discussion

Epidemiological Discussion

Our results obtained with k-support regularized logistic regression and based on a nationally representative sample of the US adult population show that the three most relevant groups of variables were represented by lifetime alcohol use disorder, as a major inclusion criteria for the study, lifetime drug use disorder, and lifetime mental health problems. This is in line with recent epidemiological studies that suggest that alcohol use disorder is in 74.9% co-occurrence with illegal drug use disorder and that comorbid alcohol use disorder is associated with a two-fold increase in the likelihood of perceived unmet need for illegal drug use disorder [START_REF] Melchior | Treatment for illegal drug use disorders: the role of comorbid mood and anxiety disorders[END_REF]. Furthermore, lifetime alcohol use disorder was shown to be strongly and significantly associated with lifetime psychiatric disorders [START_REF] Hasin | Prevalence, correlates, disability, and comorbidity of DSM-IV alcohol abuse and dependence in the United States: Results from the national epidemiologic survey on alcohol and related conditions[END_REF]. In addition, people with co-occurring substance and mental health problems tend to prefer to use mental health services instead of specialized substance use programs [START_REF] Melchior | Treatment for illegal drug use disorders: the role of comorbid mood and anxiety disorders[END_REF].

Family history of substance (alcohol, drugs) use disorders and behavioral problems was the fourth most relevant group of variables. The heritability of alcohol consumption is estimated at 35% to 40% in twin studies [START_REF] Clifford | A genetic and environmental analysis of a twin family study of alcohol use, anxiety, and depression[END_REF]. Adoptive parents with alcohol use disorder increase the probability of their adoptive children to have the same disorder [START_REF] Cotton | The familial incidence of alcoholism: A review[END_REF]. Furthermore, the aggregation of drug dependence and alcohol dependence within some families suggests common mechanisms for these disorders [START_REF] Nurnberger | A family study of alcohol dependence: Coaggregation of multiple disorders in relatives of alcohol-dependent probands[END_REF].

The type of treatment for alcohol use disorder is the fourth most relevant group of variables selected. From an epidemiological point of view, it is important to evaluate the current state of treatment in the population with alcohol use disorders in order to give recommendations for the improvement of access to care for these group of people.

Additionally, our model includes variables on treatment receipt for generalized anxiety disorder. The co-occurrence of generalized anxiety disorder (GAD) and alcohol use disorder (AUD) disorder with the rate of comorbidity ranging from 8.3% to 56.2% is well documented [START_REF] Kushner | The relation between alcohol problems and the anxiety disorders[END_REF]. GAD was also found to be a relevant factor among individuals with alcohol use disorder seeking outpatient substance abuse treatment [START_REF] Smith | Comorbidity of generalized anxiety disorder and alcohol use disorders among individuals seeking outpatient substance abuse treatment[END_REF].

Overall, variable selection performed with the help of the k-support Regularized Logistic is soundly supported by the scientific evidence. Selected variables can be used in the study that can lead to recommendation for an improvement of access to treatment among population with alcohol use disorder.

From an epidemiological point of view, variable selection using such a method is useful for identifying the most important risk factors and disregarding variables that do not add additional information in a reliable and time efficient way, therefore limiting the number of variables required for statistical analyses and improving the reliability of the final results. The selected variables may provide a suggestion for data collection for future epidemiological studies.

Statistical Discussion

Our quantitative evaluation supports the use of k-support norm regularized logistic regression over the compared methods, 1 and 2 regularized logistic regression. The use of sparsity regularization is firmly associated with significant improved empirical performance over 2 regularization and baseline methods, and the k-support norm performed better on average than 1 regularization, though a Wilcoxon signed rank test was not able to reject the null hypothesis in the comparison between k-support and 1 regularization on this data. Nevertheless, the k-support norm is to be preferred over 1 regularization even when statistically tied due to its improved handling of confounding factors and correlated signals (see Section 3.3.4). This is of particular interest when considering data such as those from NESARC in which there are a large number of variables incorporating potential confounders and comorbidities.

Conclusions

In this paper, we introduced a novel multi-class generalization of regularized logistic regression using the k-support norm. We have shown its application to the discovery of predictors of mental health service utilization on a large scale epidemiological database, NESARC. The results were both qualitatively and quantitatively analyzed, showing statistically significant improvement of k-support norm regularized logistic regression over 2 regularized logistic regression and baseline methods, and a statistical tie between k-support regularization and 1 regularization, but with better average performance and improved statistical characteristics. Qualitative analysis of the ranked variables showed a clear trend elucidating the relationship between covariates of the primary variable of interest. Discovering factors for not obtaining help with drinking is a vital public health question, with high societal and public policy impacts, and k-support regularized logistic regression has shown to be a powerful tool for analyzing this problem.

Figure 1 :

 1 Figure1: ROC curves for k-support norm regularized logistic regression, ridge regression, and the curve obtained by selecting the best performing single variable. The ROC curve for k-support norm regularized logistic regression is substantially higher than the other two curves. The curve for the single best variable is piecewise linear due to the selection of a categorical variable.

Figure 2 :

 2 Figure 2: Wilcoxon signed rank test performed pairwise between 1 , 2 and k-support regularized variants of logistic regression. Both k-support and

Figure

  Figure3: ROC curves for k-support norm regularization, 1 regularization, and 2 regularization of logistic regression. k-support and 1 regularization perform substantially better than 2 regularization. k-support regularization performs slightly better than 1 regularization on average, but without statistical significance (see Figure2). However, k-support regularization is to be favored due to its better handling of correlated variables (see Section 3.3.4).

Figure 4 :

 4 Figure4: Empirical convergence of k-support regularized logistic regression for varying values of k. 1 regularization (k = 1) corresponds to FISTA optimization (see Section 3.3.3 for details). We use the same k-support code to compute 2 regularization (k = d), although in practice one could make use of optimization methods for twice-differentiable objective functions in this special case. There is a dependence between k and the speed of optimization, due to the longer computation required per iteration for larger values of k[START_REF] Argyriou | Sparse prediction with the k-support norm[END_REF]. All computation times are feasible even on a single core for standard commercial hardware, with convergence ranging from 1-35 seconds depending on the value of k.

Table 1 :

 1 NESARC variables referred to in the article. The first column specifies the name of the variable used in the text for brevity. Alcohol Abuse and Alcoholism (NIAAA). Wave 1 NESARC data on which this study is based were collected during 2001-2002 through computer-assisted personal interviews (CAPI) in face-to-face household settings. The sample included 43,093 respondents ages 18 and older, representing the civilian, non-institutionalized adult population in the United States, including all 50 States and the District of Columbia. Military personnel living off-base and residents in non-institutionalized group quarters housing, such as boarding houses, shelters, and dormitories, were also included. All participants were interviewed at home by experienced lay interviewers who received extensive training and supervision. All procedures, including informed consent, received full ethical review and approval from the U.S. Census Bureau and U.S. Office of Management and Budget.

	Variable NESARC Code	Full Name
	A	S2CQ4A	EVER THOUGHT SHOULD SEEK HELP WITH DRINKING BUT
			DIDN'T GO
	B	ALCABDEP12DX	ALCOHOL ABUSE/DEPENDENCE IN LAST 12 MONTHS
	C	ALCABDEPP12DX ALCOHOL ABUSE/DEPENDENCE PRIOR TO THE LAST 12 MONTHS

on

Table 2 :

 2 The highest weighted variables selected from the NESARC database by k-support regularized logistic regression.INHALANT SOLVENT ABUSE DEPENDENCE PRIOR TO THE LAST 12 MONTHS 10 PANIC DISORDER WITH AGORAPHOBIA PRIOR TO THE LAST 12 MONTHS ILLNESS INDUCED 11 AGORAPHOBIA WITH NO HISTORY OF PANIC DISORDER IN LAST 12 MONTHS ILLNESS INDUCED 12 INHALANT SOLVENT ABUSE DEPENDENCE IN LAST 12 MONTHS 13 HEROIN ABUSE DEPENDENCE IN LAST 12 MONTHS 14 AGORAPHOBIA WITH NO HISTORY OF PANIC DISORDER PRIOR TO THE LAST 12 MONTHS 15 COCAINE ABUSE DEPENDENCE IN LAST 12 MONTHS 16 TRANQUILIZER ABUSE DEPENDENCE IN LAST 12 MONTHS 17 ADOPTIVE MOTHER EVER AN ALCOHOLIC OR PROBLEM DRINKER 18 AMPHETAMINE ABUSE DEPENDENCE PRIOR TO THE LAST 12 MONTHS 19 AMPHETAMINE ABUSE DEPENDENCE IN LAST 12 MONTHS 20 HALLUCINOGEN ABUSE DEPENDENCE IN LAST 12 MONTHS 21 HALLUCINOGEN ABUSE DEPENDENCE PRIOR TO THE LAST 12 MONTHS 22 PANIC DISORDER WITH AGORAPHOBIA IN LAST 12 MONTHS ILLNESS INDUCED AND 23 OPIOID ABUSE DEPENDENCE IN LAST 12 MONTHS 24 MANIC DISORDER IN LAST 12 MONTHS ILLNESS INDUCED AND SUBSTANCE INDUCED 25 COCAINE ABUSE DEPENDENCE PRIOR TO THE LAST 12 MONTHS 26 WENT TO COUNSELOR THERAPIST DOCTOR OTHER PERSON FOR HELP TO IMPROVE MOOD 27 STAYED OVERNIGHT IN HOSPITAL BECAUSE OF DYSTHYMIA 28 WENT TO EMERGENCY ROOM FOR HELP BECAUSE OF DYSTHYMIA 29 DOCTOR PRESCRIBED MEDICINE DRUG TO IMPROVE MOOD MAKE YOU FEEL BETTER 30 DYSTHYMIA IN LAST 12 MONTHS ILLNESS INDUCED AND SUBSTANCE INDUCED RULED OUT 31 MANIC DISORDER PRIOR TO THE LAST 12 MONTHS ILLNESS INDUCED AND 32 GENERALIZED ANXIETY IN LAST 12 MONTHS ILLNESS INDUCED AND SUBSTANCE INDUCED 33 BLOOD NATURAL MOTHER EVER AN ALCOHOLIC OR PROBLEM DRINKER 34 CANNABIS ABUSE DEPENDENCE IN LAST 12 MONTHS 35 BLOOD NATURAL MOTHER EVER HAD PROBLEMS WITH DRUGS 36 ALCOHOL ABUSE DEPENDENCE PRIOR TO THE LAST 12 MONTHS 37 DYSTHYMIA PRIOR TO THE LAST 12 MONTHS ILLNESS INDUCED AND SUBSTANCE INDUCED 38 HYPOMANIC DISORDER IN LAST 12 MONTHS ILLNESS INDUCED AND SUBSTANCE INDUCED 39 ANY FULL BROTHERS EVER ALCOHOLICS OR PROBLEM DRINKERS 40 WENT TO COUNSELOR THERAPIST DOCTOR OTHER PERSON FOR HELP WITH PANIC ATTACKS 41 WENT TO EMERGENCY ROOM FOR HELP BECAUSE OF PANIC ATTACKS 42 STAYED OVERNIGHT IN HOSPITAL BECAUSE OF PANIC ATTACKS 43 DOCTOR PRESCRIBE MEDICINE DRUG FOR PANIC ATTACKS 44 NATURAL GRANDMOTHER ON MOTHER S SIDE EVER AN ALCOHOLIC OR PROBLEM DRINKER 45 DOCTOR EVER PRESCRIBE MEDICINE DRUG FOR GENERALIZED ANXIETY

	NESARC variable name
	1 EVER SOUGHT HELP BECAUSE OF DRINKING
	2 ANY FULL SISTERS EVER ALCOHOLICS OR PROBLEM DRINKERS
	3 SEDATIVE ABUSE DEPENDENCE PRIOR TO THE LAST 12 MONTHS
	4 TRANQUILIZER ABUSE DEPENDENCE PRIOR TO THE LAST 12 MONTHS
	5 OPIOID ABUSE DEPENDENCE PRIOR TO LAST THE 12 MONTHS
	6 HEROIN ABUSE DEPENDENCE PRIOR TO THE LAST 12 MONTHS
	7 SEDATIVE ABUSE DEPENDENCE IN LAST 12 MONTHS
	8 ADOPTIVE FATHER EVER AN ALCOHOLIC OR PROBLEM DRINKER
	9

regularization perform significantly better than

regularization. k-support regularization has a higher AUC than 1 regularization on average.

Acknowledgements

This work was supported in part by the European Commission through ERC Grant 259112, and FP7-MC-CIG 334380. The funding agencies had no involvement in the decision to submit the article for publication.


V Machairas

E Decencière

T Walter

WATERPIXELS: SUPERPIXELS BASED ON THE WATERSHED TRANSFORMATION

Keywords: Superpixels, watershed, segmentation

HAL is

INTRODUCTION

"Superpixels" (SP) are homogeneous regions resulting from a low-level segmentation of an image and typically acting as primitives for further analysis such as detection, segmentation and classification of objects. Superpixels should have the following properties:

1. homogeneity: pixels of a given SP should present similar colors or gray levels;

2. connected partition: each SP is made of a single connected component; SP do not overlap;

3. adherence to object boundaries: object boundaries should be included in SP boundaries;

4. regularity: SPs should form a regular pattern on the image. This property is often desirable as it makes the SP more convenient to use for subsequent analysis steps.

It is clear, that in practice, properties [START_REF] Ren | Learning a classification model for segmentation[END_REF] and ( 4) are somewhat contradictory and consequently a good solution typically aims at finding a compromise between these two requirements.

Low-level segmentations have been used for a long time as first step towards segmentation [START_REF] Monga | An optimal region growing algorithm for image segmentation[END_REF][START_REF] Marcotegui | Bottom-up segmentation of image sequences for coding[END_REF]. The term superpixel was coined much later [START_REF] Ren | Learning a classification model for segmentation[END_REF] in this context, albeit in a more constrained framework. This approach has raised increasing interest since then.

Various methods exist to compute SP, most of them based on graphs, geometrical flows or k-means (see [START_REF] Levinshtein | TurboPixels: Fast Superpixels Using Geometric Flows[END_REF], [START_REF] Veksler | Superpixels and supervoxels in an energy optimization framework[END_REF], [START_REF] Zeng | Structure-sensitive superpixels via geodesic distance[END_REF] and [START_REF] Achanta | SLIC superpixels compared to state-ofthe-art superpixel methods[END_REF]).

In principle, the watershed transformation is well suited for SP generation: first because it gives a good adherence to object boundaries when computed on the gradient of the image, second because it allows to control the number and spatial arrangement of the resulting regions by the choice of markers and third because it offers linear complexity with the number of pixels in the image. Indeed, it has been used to produce low-level segmentations in several applications, including computation intensive 3D applications [START_REF] Andres | Segmentation of SBFSEM Volume Data of Neural Tissue by Hierarchical Classification[END_REF][START_REF] Stawiaski | Interactive liver tumor segmentation using watershed and graph cuts[END_REF], in particular when shape regularity of the elementary regions was not required.

Previous publications claimed that the watershed transformation does not allow for the generation of spatially regular SP [START_REF] Levinshtein | TurboPixels: Fast Superpixels Using Geometric Flows[END_REF][START_REF] Achanta | SLIC superpixels compared to state-ofthe-art superpixel methods[END_REF]. Here, we show that on the contrary, a correctly used watershed can lead to efficient SP computation. We name this method of generating SP "waterpixels".

WATERPIXELS GENERATION METHOD

Let f : D → V be an image, where D is a bounded connected subset of Z 2 , and V a set of values, typically {0, . . . , 255} when f is a grey level image, or {0, . . . , 255} 3 for colour images.

As any watershed based segmentation, waterpixels are based on two principles: the definition of markers, from which the flooding starts, and the definition of a gradient (the image to be flooded), as illustrated in Figure 2.

The markers and the gradient

First, we choose a set of N points {v i } 1≤i≤N in D, called cell centers. As we want some degree of spatial regularity, they can be placed on the vertices of a square or hexagonal grid. Given a distance on D, we denote by σ the distance between closest points, which is linked to N . It is worth noting that the distance does not have to be Euclidean.

A Voronoi tesselation allows to associate to each v i a Voronoi cell (see Fig. 2.c). For each such cell, a homothety centered on v i with factor ρ (0 < ρ ≤ 1) leads to the computation of the final cell C i . This last step allows for the creation of a margin between neighbouring cells, in order to avoid the selection of minima too close from each other.

As each cell is meant to correspond to the generation of a unique waterpixel, our method, through the choice of cell centers, offers total control over the number of SP, with a strong impact on their size and shape if desired.

Second, a gradient image g is computed from f (see respectively Fig. 2.a and Fig. 2.b as an example). The choice of the gradient operator depends on the image type, e.g. for grey level images we might choose a morphological gradient. Within each cell, a single minimun of g will be used as marker. If several minima are present, then the one with the highest volume extinction value [START_REF] Vachier | Extinction values: A new measurement of persistence[END_REF] is used (see Fig. 2

.d).

If no minimum is present, the center of the cell is used as marker. Final selection of the markers is illustrated in Fig. 2.g.

Spatial regularization of the gradient

When directly flooding the gradient image g from the selected minima, the resulting region often suffers from irregular borders. Here, we propose to control this irregularity by using a spatially regularized gradient (see Fig. 2.e):

g reg = g + k 2d σ ( 1
)
where d is the above introduced distance function to the cell centers, i.e. d(p) is the distance of pixel p to the closest cell center, and k is the spatial regularization parameter, which takes its values within + . Note that d is normalized with the grid radius, in order to make k independent from the choice of σ. Resulting waterpixels are presented in Fig. 2.f. The choice of k is application dependent: when k equals zero, no regularization of the gradient is applied; when k → ∞, we approach the regular grid. This behaviour is illustrated in Fig. 2.h, 2.i and 2.j, for k equal to 0, 4 and 10 respectively.

APPLICATION TO THE BERKELEY SEGMENTATION DATABASE

In order to evaluate the proposed strategy, we have applied it to the Berkeley segmentation database [START_REF] Martin | A database of human segmented natural images and its application to evaluating segmentation algorithms and measuring ecological statistics[END_REF], and compared it with a state-of-the art method.

Implementation

We have found (data not shown) that it is beneficial to preprocess the images from the database using an area opening followed by an area closing, both of size σ 2 /16. This operation efficiently removes details which are clearly smaller than the expected waterpixel area and which should therefore not give rise to a superpixel contour. The Lab-gradient is adopted here in order to best reflect our visual perception of color differences and hence the pertinence of detected objects. The cell centers correspond to the vertices of a square regular grid of step σ. A square grid has been chosen to make comparison with SLIC (see next section), which is also based on a square grid, easier. For the distance map, we have used the distance metric d(x, y) = x -y ∞ which corresponds to the placement of markers on a rectangular grid. The margin parameter ρ is set to 2/3. Note that the gradient images have integer values, which enables to use a fast implementation of the watershed transformation based on a hierarchical queue.

Evaluation criteria

Superpixel methods produce an image partition. In order to compute the superpixels borders, we use a morphological gradient. Note that the resulting contours are two pixels wide. To this set, we add the image contours. The final set is denoted C. A ground truth image GT corresponds to contours of the objects to be segmented.

In order to quantitatively assess the quality of our superpixels, we have used two evaluation criteria:

• Boundary-recall (BR), which measures adherence to boundaries without penalizing over-segmentation, and is defined as the percentage of ground-truth contour pixels GT which fall within less than 3 pixels from superpixel boundaries: • Contour density (CD), which measures the irregularity of the partition, and is defined as the ratio between the number of contour pixels of SP and the total number |D| of pixels in the image:

BR = |{p ∈ GT, d(p, C) < 3}| |GT | (2)
CD = |C| |D| (3)
The proposed method is benchmarked against the stateof-the-art method simple linear iterative clustering. SLIC [START_REF] Achanta | SLIC superpixels compared to state-ofthe-art superpixel methods[END_REF] is a method based on k-means which offers a linear complexity with respect to the number of pixels in the image. The properties of the SP can be controlled by two parameters: n (number of SP) and m which influences the distance metric used and thereby the spatial regularity. SLIC is currently considered as one of the best methods to compute regular SP in terms of quality as well as computation time, but has a major drawback as it needs post-processing to treat disconnected superpixels.

Both methods have been applied on the subset "train" of the Berkeley segmentation database, containing 300 images of sizes 321 × 481 or 481 × 321 pixels. Approximately 6 human-annotated ground-truth segmentations are given for each image. These ground-truth images correspond to manually drawn contours.

Benchmark

Results for boundary-recall and contour density, expressed as a function of the number of superpixels in the image, and averaged over the whole database, are shown in figure 3. Blue and red curves correspond to varying regularization parameters k and m respectively for waterpixels and SLIC.

From figure 3 (a) and (b) we see that the working points tested for the two algorithms were different in most cases. In order to make a fair comparison, we choose a pair of parameters (m, k) such that they show similar boundary recall and compare their behavior for contour density, e.g. m = 20, k = 4. We see that the results for contour density are indeed very similar. We therefore conclude that both methods achieve comparable results.

DISCUSSION

While results are comparable, there is one major difference in the construction of the algorithm: the SLIC approach does not impose any connectivity constraint. The resulting superpixels are therefore not necessarily connected, which requires some ad hoc postprocessing step. In contrast, waterpixels are connected by definition, and the connectivity constraint is actually implemented in the distance used.

Furthermore, waterpixels offer a nice perspective to efficiently build hierarchical partitions based on superpixels. Indeed, the computation of the watershed naturally produces a The same graph can be used to build waterpixels at different resolutions, or to fuse waterpixels in order to obtain a high-level segmentation.

CONCLUSION AND PERSPECTIVES

We have proposed a general strategy to build superpixels based on the watershed segmentation method. We have shown that it performs as well as state-of-the-art methods on the Berkeley segmentation database. Its theoretical complexity is also linear with respect to the number of image pixels, as for SLIC.

The current implementation of our method, as other superpixel methods, depends on two main parameters σ and k. However, two "hidden" parameters, linked to minima selection, need to be tuned: the cells margins and the preprocessing area parameter. On the other hand, we have not yet fully exploited the main degree of freedom the method proposes which consists in the placement of markers. We are therefore currently investigating the possibility to place markers optimally, such that the preprocessing steps might become unnecessary.

Fig. 2 .

 2 Fig. 2. Illustration of waterpixels generation: (a): gradient of the original image (b); (c): regular grid with square cells and step σ=50 pixels; (d) process for selecting one minimum of (a) per cell of (c) -final markers set is presented in (g); (e): regularized gradient (k=10); (d) Resulting waterpixels obtained by applying the watershed transformation on (e) with markers of (g). Impact of the spatial regularization on the regularity of resulting waterpixels : (h) k=0, (i) k=4 and (j) k=10.

Fig. 3 .

 3 Fig. 3. Benchmark: performance comparison of Waterpixels to SLIC