Analysis of overfitting in the regularized Cox model

1 Introduction

Inference of parameters for generalized linear models using the maximum likelihood (ML) protocol becomes increasingly biased due to overfitting as the ratio $ζ = p / N$ increases, where $p$ is the number of covariates and $N$ the number of training data. Overfitting occurs when model parameters seek to explain not only the ‘signal’ but also the ‘noise’ in training data, and is characterized by a difference in outcome prediction accuracy between training and validation samples. See e.g. [1, 3, 4, 5]

for examples from logistic regression, gamma distributions and Cox models

[1]. Hence, standard statistical significance tests for regression coefficients, being usually based on asymptotic results derived for fixed

p

, become increasingly inaccurate [6]. Unfortunately, in post-genome medicine, having large ratios

ζ

is the rule rather than the exception. This prompted epidemiologists to formulate heuristic rules for avoiding overfitting, such as limits on the number of events per variable

[7, 8, 9, 10, 11]. Approximate recipes for correcting ML estimates were developed in e.g. [12, 13]

. Alternative methods of addressing the overfitting problem include feature selection and regularization. In feature selection one seeks to identify a subset of covariates that are informative of outcomes

[14, 15, 16]

. Its advantages include reduction in the required computational resources, and increased interpretability. In regularization one adds a penalty term to the objective function of ML inference (which can alternatively be derived from a prior in Bayesian inference) to suppress the number or magnitude of the model parameters

[17, 18]. Application of regularization to survival analysis with high-dimensional covariates is studied widely, see e.g. [19, 20] and references therein.

Figure 1: Synthetic data with dimension $p = 25$ and $N = {10, 100, 1000}$ are generated using the logistic regression model. The ML estimate of model parameters is found numerically using the Nelder-Mead algorithm. The overfitting measure $E$ is plotted after each iteration. The starting value model parameters $β$
in the minimization search is the zero vector, giving a positive value of E (implying an underfitted model).

Figure 2: Comparison of true and inferred regression coefficients for the Cox proportional hazards model. A systematic bias is found for non-zero values of $ζ$ , which can be corrected with regularization. Synthetic survival data were generated [21] using Gaussian covariates ( $p = 500, N = 833, ζ = 0.6$ ). The regression coefficients are inferred [17]
using ML (no regularization) or Maximum A Posteriori Probability regression (MAP, with regularization). Data points on the diagonal imply perfect inference.

A recent study [1] provided a new approach to overfitting in survival analysis. It showed how the replica method from statistical physics can be used to model ML inference analytically as the zero noise limit of a suitably defined stochastic minimization, starting from an information-theoretic measure of overfitting. The theory predicted the quantitative relation between ML-inferred and true parameters in the Cox model [1]

, and a phase transition at

ζ = 1

Let us denote the set of model parameters as $ϑ$ , and the data as $D$

. The observation that ML inference is equivalent to minimization of the Kullback-Leibler divergence between the empirical data distribution

{^P}_{D}

and the parametrized distribution

P_{ϑ}

assumed as a model of the data, suggests [1] using

E (ϑ, D) \equiv D ({^P}_{D} ∥ P_{ϑ}) - D ({^P}_{D} ∥ P_{ϑ^{⋆}})

as a measure of overfitting¹¹1

A similar idea for comparing estimators of probability distributions was used in

[22], using the Lévy distance rather than the KL divergence. Other measures of overfitting can be found in e.g. [20, 23]., in which

ϑ^{⋆}

are the true (but a priori unknown) parameter values. Perfect regression implies

E = 0

, underfitting implies

E > 0

, and overfitting implies

E < 0

. To gain more intuition for this measure, we generate synthetic data from a simple logistic regression model, find the ML estimators of its parameters, and calculate

E

. Here the parameters are

{β_{μ}}_{μ = 0}^{p}

, the data are

D = {(x_{1}, y_{1}), \dots, (x_{N}, y_{N})}

, with

x_{i} \in R^{p + 1}

and

y_{i} \in {0, 1}

, and we use the short-hand

β \cdot x = \sum_{μ = 0}^{p} β_{μ} x_{μ}

(with the convention

x_{0} = 1

). The outcome likelihood

P_{β} (y | x)

and the measure

E

are given by

	$P_{β} (y \| x)$	$=$	$(11+\rme−β⋅x)y(11+\rmeβ⋅x)1−y$		(1.1)
	$E (β^{⋆}, D)$	$=$	$1NN∑i=1{yilog(1+\rme−β⋅xi1+\rme−β⋆⋅xi)+(1−yi)log(1+\rmeβ⋅xi1+\rmeβ⋆⋅xi)}$		(1.2)

Results are shown in Figure 1. When $ζ = 0.025$ , $E$ converges towards zero during the minimization, indicating perfect parameter recovery. As the number of samples in the data set is reduced, giving $ζ = 0.25$ and $ζ = 2.5$ , $E$ converges to increasingly negative values. Since there is no model mismatch (the data were generated from a logistic model), the negative values of $E$ indicate overfitting.

Switching from maximum likelihood to maximum a posteriori estimators implies adding a penalty term to the likelihood:

D ({^P}_{D} ∥ P_{ϑ}) \to D ({^P}_{D} ∥ P_{ϑ}) - log p (ϑ)

where

p (ϑ)

represents a parameter prior, giving

	$E (ϑ^{⋆}, D)$	$\equiv$			(1.3)
		$=$	$min ϑ {\frac{1}{N} N \sum i = 1 log \frac{p (t_{i} \| z_{i}, ϑ^{⋆}) p (ϑ^{⋆})}{p (t_{i} \| z_{i}, ϑ) p (ϑ)}}$		(1.3)

MAP regression is equivalent to minimizing the quantity (1.3). This minimization should in principle be over all $ϑ$ , but may in practice be constrained to simplify the calculation (see e.g. [24, 25, 26]). For generalized linear models, commonly used priors are $p (β) \propto exp [- η \sum_{μ = 1}^{p} | β_{μ} |]$ (giving $L_{1}$ regularization²²2This choice promotes sparsity in the regression coefficient vector $β$ , which would result in a horizontal line segment passing through the origin in Fig. 2. Since our theory aims to predict the slope of the data clouds in Fig. 2, we will not pursue $L_{1}$ regularizers in this paper.

, or ‘LASSO’ regression

[18]) and

p (β) \propto exp [- η \sum_{μ = 1}^{p} β_{μ}^{2}]

(giving

L_{2}

regularization, or ‘ridge’ regression).

In the present paper we generalize the replica analysis of [1] from ML to MAP inference, upon adding an $L_{2}$ regularization term to the log-likelihood function. This term suppresses overfitting effects, and removes the ML phase transition of the Cox model [1] at $ζ = 1$ ; see e.g. Fig. 2. In the presence of an $L_{2}$ regularizer, correlations between covariates can no longer be transformed away, as was done in [1], leading to the appearance in the theory of the population covariance matrix $A$

of the covariates. Under mild restrictions on the eigenvalue spectrum of this matrix, we show how an accurate theory of overfitting for the regularized Cox proportional hazards model can be developed, in spite of such additional mathematical complications, including for the previously inaccessible regime

ζ > 1

. We find, as in [1], that the replica symmetric version of the theory is sufficient to explain accurately the behaviour of interest. The resulting equations can also be used to predict the amount of regularization needed for unbiased regression, expressed in term of spectum of

A

and the ratio

ζ

2 Replica analysis of regularized Cox regression

2.1 Generalized replica formalism to include priors

Following [1], we interpret minimization of (1.3

) as computing the ground state energy of a statistical mechanical system with degrees of freedom

ϑ

and Hamiltonian

H (ϑ | ϑ^{⋆}, D)

, at inverse temperature

γ

, where

D = {(t_{1}, z_{1}), \dots, (t_{N}, z_{N})}

and

H (ϑ | ϑ^{⋆}, D)

=

log N \prod i = 1 [\frac{p (t_{i} | z_{i}, ϑ^{⋆}) p (ϑ^{⋆})}{p (t_{i} | z_{i}, ϑ) p (ϑ)}]

(2.1)

We define the associated free energy, which we average over the disorder (the microscopic realization of $D$ ), and can compute the disorder-averaged ground state energy as the $γ \to \infty$ limit of the disorder-averaged energy density $E_{γ} (ϑ^{⋆})$ , where

	$E_{γ} (ϑ^{⋆})$	$=$	$−1N∂∂γ⟨log∫\rmdϑ N∏i=1[p(ti\|zi,ϑ)p(ϑ)p(ti\|zi,ϑ⋆)p(ϑ⋆)]γ⟩D$

The replica identity $⟨ log Z ⟩ = {lim}_{n \to 0} n^{- 1} log ⟨ Z^{n} ⟩$ is subsequently used to simplify the average of the logarithm (see [1] for details and references), giving in the present case

$E_{γ} (ϑ^{⋆})$	$=$	$−∂∂γlimn→01Nnlog⟨{∫\rmdϑ N∏i=1[p(ti\|zi,ϑ)p(ϑ)p(ti\|zi,ϑ⋆)p(ϑ⋆)]γ}n⟩D$	(2.3)
	$=$
	$=$
		$×{∫\rmdz\rmdt p(z)p(t\|z,ϑ⋆)n∏α=1[p(t\|z,ϑα))p(t\|z,ϑ⋆)]γ}N$

Equation (2.3

) is applicable to any parametric model

p (t | z, ϑ)

and any prior

p (ϑ)

. See also [27] for alternative results on the use of the replica method in statistical inference. We will now make a specific choice for

p (t | z, ϑ)

, and use (2.3) to develop a theory for regression and overfitting in regularized Cox models with Gaussian priors.

2.2 Application to the regularized Cox proportional hazards model

Cox’s proportional hazards model originally described in [2] assumes a parametrization of the form

p (t | z, ϑ)

=

λ(t)\rmeβ⋅z−exp(β⋅z)∫t0\rmdt′ λ(t′)

(2.4)

Its parameters are the coefficients $β \in R^{p}$ , and a base hazard rate $λ (t)$ (a nonnegative function defined for $0 \leq t < \infty$ ). Substituting $ϑ = {β, λ}$ translates (2.3) into

	$E_{γ} (β^{⋆}, λ^{⋆})$	$=$	$−∂∂γlimn→01Nnlog∫{\rmdλ1…\rmdλn}∫\rmdβ1…\rmdβn{n∏α=1[p(βα)p(β⋆)]γ}$		(2.5)
			$×{∫\rmdz\rmdt p(z)p(t\|z,β⋆,λ⋆)n∏α=1[p(t\|z,βα,λα))p(t\|z,β⋆,λ⋆)]γ}N$		(2.5)

Functional integrals are written as $∫{\rmdλ}$ , the true parameters responsible for the data are written as ${β^{⋆}, λ^{⋆}}$ , and we follow the standard convention for regularized Cox models of only including a prior for the association parameters (equivalently, assuming an improper, or ‘flat’, prior for the base hazard rate). Our $L_{2}$ prior is $p (β) \propto exp (- η β^{2})$ . To proceed with the analytical treatment, we assume that the covariate vectors $z_{i}$ are drawn independently from a population distribution with zero mean and covariance matrix $A$ . The introduction of regularization means that the regression equations for correlated covariates can no longer be transformed to those corresponding to uncorrelated ones. This leads to a more complex theory than [1], and ultimately to conditions on the eigenvalue spectrum of $A$ .

Our analysis is carried out in the regime where both $N, p \to \infty$ but with fixed ratio $ζ = p / N \sim O (1)$ . To retain non-zero event times, even for $p \to \infty$ , we must rescale the regression coefficients according to $β \cdot z \to β \cdot z / \sqrt{p}$ , resulting in $β \cdot z \sim O (1)$ . Without this rescaling we would have event time distributions with all weight concentrated on $t \to 0$ and $t \to \infty$ . We also replace $β^{⋆}$ by $β^{0}$ , to allow for more compact notation. Following [1] we next introduce

p(y|β0,…,βn)=∫\rmdz p(z)n∏α=0δ[yα−βα⋅z√p]

(2.6)

where $y = {y^{0}, y^{1}, \dots, y^{n}} \in R^{n + 1}$ . The magnitude of $\rmeβ⋅z$ represents the relative risk of failure, compared to that of an ‘average’ individual (with $z = 0$ ). Therefore $y$ can be considered a vector of risk scores. Our energy density then becomes

	$E_{γ} (β^{⋆}, λ^{⋆})$	$=$	$−∂∂γlimn→01Nnlog∫{\rmdλ1…\rmdλn}∫\rmdβ1…\rmdβnn∏α=1[p(βα)p(β0)]γ$		(2.7)
			$×{∫\rmdy p(y\|β0,…,βn)∫\rmdt p(t\|y0,λ0)n∏α=1[p(t\|yα,λα)p(t\|y0,λ0)]γ}N$		(2.7)

in which now

p (t | y, λ)

=

λ(t)\rmey−exp(y)∫t0\rmdt′ λ(t′)

(2.8)

To proceed we assume that $p (y | β^{0}, \dots, β^{n})$ is Gaussian. This holds for any $N$ and $p$ as soon as $p (z)$

is Gaussian, and for non-Gaussian covariate statistics it will generally hold due to the Central Limit Theorem if the correlations among the covariates are weak, and

N

and

p

are large. Since we assumed

∫\rmdz p(z)z=0

, the risk score distribution is now given by

p(y|β0,…,βn)=\rme−12y⋅C−1[{β}]y√(2π)n+1detC[{β}]

(2.9)

It is determined in full by the $(n + 1) \times (n + 1)$ covariance matrix $C [{β}]$ , with entries

C_{α ρ} [{β}]

=

∫\rmdz p(z)(βα⋅z√p)(βρ⋅z√p) =1pβα⋅Aβρ

(2.10)

The entries of $A$ are given by $Aμν=∫\rmdz p(z)zμzν$ . The ${C_{α ρ} [{β}]}$ measure the similarity between the $p$ -dimensional vectors formed by the regression parameters in different replicas. For each replica pair $(α, ρ)$ we use the integral representation of the Dirac delta function, and rescale the conjugate integration parameter by $p$ ,

1=∫\rmdCαρ δ[Cαρ−1pβα⋅Aβρ]=∫\rmdCαρ\rmd^Cαρ2π/p\rme\rmip^Cαρ(Cαρ−βα⋅Aβρ)

(2.11)

in order to simplify expression (2.7) to

$E_{γ} (β^{⋆}, λ^{⋆})$	$=$	$−∂∂γlimn→01Nnlog∫{\rmdλ1…\rmdλn}∫\rmdC\rmd^C \rme\rmip∑nα,ρ=0^CαρCαρ(2π/p)(n+1)2$	(2.12)

		$×∫\rmdβ1…dβn \rme−ηγ∑nα=1[(βα)2−(β0)2]−\rmi∑nα,ρ=0^Cαρβα⋅Aβρ$

The quadratic nature of the exponent in the $β$ integral, a consequence of having chosen $L_{2}$ regularization, allows for a closed form solution. Changing the penalty term to $L_{1}$ or $L_{q}$ with $q > 2$ would significantly complicate the integrals.

2.3 Conversion into a saddle point problem

With a modest amount of foresight we transform $^C=−12\rmiD$ , and introduce the short-hand $~ β \equiv A^{\frac{1}{2}} β$ . To evaluate the Gaussian $β$ integral in (2.12) we define the $n p \times n p$ matrix $Ξ$ and the $n p$ -dimensional vector $ξ$ , with entries

Ξ_{α μ; β ν} = 2 η γ δ_{α β} (A^{- 1})_{μ ν} + δ_{μ ν} D_{α β}, ξ_{μ}^{α} = - D_{0 α} {~ β}_{μ}^{0}

(2.13)

With these definitions we may write the Gaussian integral in (2.12) as

			$∫(n∏α=1\rmd~βα\rme−ηγ~βα⋅A−1~βα)\rme−12∑nα,ρ=1Dαρ~βα⋅~βρ−∑nρ=1D0ρ~β0⋅~βρ$		(2.14)
		$=$	$\rme12ξ⋅Ξ−1ξ∫\rmd~β \rme−12(~β−Ξ−1ξ)⋅Ξ(~β−Ξ−1ξ)=(2π)np2√detΞ\rme12ξ⋅Ξ−1ξ$		(2.14)

Let ${a_{μ}}$ and ${b_{α}}$ denote the eigenvalues of $A$ and $D$ , respectively. The two terms $P$ and $Q$ of the matrix $Ξ$ , with components $P_{α μ, β ν} = 2 η γ δ_{α β} (A^{- 1})_{μ ν}$ and $Q_{α μ, β ν} = δ_{μ ν} D_{α β}$

, clearly commute. The complete set of eigenvectors of

Ξ

can therefore be written as

{^uμα}

, with components

{^u}_{ν ρ}^{μ α} = u_{ρ}^{α} v_{ν}^{μ}

, and where

\sum_{ρ \leq n} D_{λ ρ} u_{ρ}^{α} = b_{α} u_{ρ}^{λ}

and

\sum_{ν \leq p} A_{λ ν} v_{ν}^{μ} = a_{μ} v_{λ}^{μ}

, and where both are normalised according to

\sum_{ρ \leq n} (u_{ρ}^{α})^{2} = \sum_{ν \leq p} (v_{ν}^{μ})^{2} = 1

. The eigenvalues of

Ξ

are then

ξ_{μ α} = 2 η γ / a_{μ} + b_{α}

, and

det Ξ = p \prod μ = 1 n \prod α = 1 (\frac{2 η γ}{a_{μ}} + b_{α}), (Ξ^{- 1})_{α μ, α^{'} μ^{'}} = n \sum β = 1 p \sum ν = 1 \frac{u_{α}^{β} v_{μ}^{ν} u_{α^{'}}^{β} v_{μ^{'}}^{ν}}{2 η γ / a_{ν} + b_{β}}

(2.15)

Hence the integral (2.14) can be written as

(2π)np2√detΞ \rme12ξ⋅Ξ−1ξ

=

\rme12nplog(2π)−12np⟨log(2ηγ/a+b)⟩+12np⟨(ξ⋅^u)2(2ηγ/a+b)−1⟩

(2.16)

where the averages in the exponents are over the eigenvalues and orthonormal eigenvectors of $Ξ$ , i.e. $⟨ f (a, b,^u) ⟩ = (n p)^{- 1} \sum_{μ = 1}^{p} \sum_{α = 1}^{n} f (a_{μ}, b_{α}, {^u}^{μ α})$ . Since $p = ζ N$ with $ζ > 0$ , the integrals over $C$ , $^C$ and the base hazard rates in (2.12) can for $N \to \infty$ be evaluated by steepest descent, provided the limits $n \to 0$ and $N \to \infty$ commute. Expression (2.16) then enables us to write the result as

lim N \to \infty E_{γ} (β^{⋆}, λ^{⋆})

=

\frac{\partial}{\partial γ} lim n \to 0 \frac{1}{n} extr Ψ (C, D, λ^{1} \dots λ^{n})

(2.17)

in which

$Ψ (C, D, λ^{1} \dots λ^{n})$	$=$	$- \frac{1}{2} ζ [n \sum α, ρ = 0 D_{α ρ} C_{α ρ} - \frac{1}{p} D_{00} ({~ β}^{0})^{2}] + \frac{1}{2} (n + 1 - n ζ) log (2 π)$	(2.18)

		$−log∫\rmdy \rme−12y⋅C−1y∫\rmdt p(t\|y0,λ0)n∏α=1[p(t\|yα,λα)p(t\|y0,λ0)]γ$

where we have defined $S^{2} = {lim}_{p \to \infty} p^{- 1} (β^{0})^{2}$ . Differentiating $Ψ$ with respect to $D_{00}$ removes $D_{00}$ from the problem, and gives $C_{00} = p^{- 1} β^{0} \cdot A β^{0} \equiv {~ S}^{2}$ .

3 Replica symmetric theory

3.1 Replica symmetric saddle points

To proceed, we make the replica symmetric ansatz, which implies assuming ergodicity of the stochastic regression process, and translates into invariance of all order parameters under all permutations of the replicas ${1, \dots, n}$ . Now, for all $1 \leq α, ρ \leq n$ :

λ^{α} (t) = λ (t), \begin{matrix} C_{0 α} = c_{0} D_{0 α} = d_{0} \end{matrix}, \begin{matrix} C_{α ρ} & = C δ_{α ρ} + c (1 - δ_{α ρ}) D_{α ρ} & = D δ_{α ρ} + d (1 - δ_{α ρ}) \end{matrix}

(3.1)

Both $C$ and $D$ are positive definite, so $C > c$ and $D > d$ . We may now write

C = ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ \begin{matrix} C_{00} & c_{0} & \dots & \dots & c_{0} c_{0} & C & c & \dots & c ⋮ & c & C & \dots & c ⋮ & ⋮ & ⋮ & ⋱ & ⋮ c_{0} & c & c & \dots & C \end{matrix} ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠, C^{- 1} = ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ \begin{matrix} B_{00} & b_{0} & \dots & \dots & b_{0} b_{0} & B & b & \dots & b ⋮ & b & B & \dots & b ⋮ & ⋮ & ⋮ & ⋱ & ⋮ b_{0} & b & b & \dots & B \end{matrix} ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠

The eigenvalues and eigenvectors of $C$ , $C^{- 1}$ and $D$ are found in [1]. $C$ has two nondegenerate eigenvalues $λ_{\pm}$ with $λ_{+} λ_{-} = [C + (n - 1) c] C_{00} - n c_{0}^{2}$ , and a further $n - 1$ fold degenerate eigenvalue $λ_{0} = C - c$ . Hence

$log det C$	$=$	$log ([C + (n - 1) c] C_{00} - n c_{0}^{2}) + (n - 1) log (C - c)$
	$=$	$log C_{00} + n log (C - c) + \frac{n (c - c_{0}^{2} / C_{00})}{C - c} + O (n^{2})$

The entries of $C^{- 1}$ are found to be

	$B_{00} = \frac{C + (n - 1) c}{C_{00} [C + (n - 1) c] - n c_{0}^{2}}, b_{0} = - \frac{c_{0}}{C_{00} [C + (n - 1) c] - n c_{0}^{2}}$		(3.3)
	$B = b + \frac{1}{C - c}, b = \frac{c_{0}^{2} - c C_{00}}{(C_{00} [C + (n - 1) c] - n c_{0}^{2}) (C - c)}$		(3.4)

Hence

	$y \cdot C^{- 1} y$	$=$	$B_{00} (y^{0})^{2} + (B - b) n \sum α = 1 (y^{α})^{2} + b (n \sum α = 1 y^{α})^{2} + 2 b_{0} y^{0} n \sum α = 1 y^{α}$

Next we turn to terms in (2.18) that involve the spectrum of $D$ . This matrix has one eigenvalue $D + (n - 1) d$ with eigenvector $v = (1, \dots, 1)$ , and the $n - 1$ fold degenerate eigenvalue $D - d$

with eigenspace

\sum_{i = 1}^{n} v_{i} = 0

. Hence

	$⟨ log (\frac{2 η γ}{a} + b) ⟩$	$=$	$\frac{1}{n p} p \sum μ = 1 [(n - 1) log (\frac{2 η γ}{a} + D - d) + log (\frac{2 η γ}{a} + D - d + n d)]$		(3.6)
		$=$	$⟨ log (\frac{2 η γ}{a} + D - d) ⟩ + ⟨ \frac{d a}{2 η γ + (D - d) a} ⟩ + O (n)$		(3.6)

Similarly, using the RS form of $ξ_{μ}^{α} = - d_{0} (A^{\frac{1}{2}} β^{0})_{μ}$ , we may write

$⟨ \frac{(ξ \cdot^u)^{2}}{2 η γ / a + b} ⟩$	$=$	$\frac{1}{n p} p \sum μ = 1 (p \sum ν = 1 n \sum ρ = 1 (A^{\frac{1}{2}} β^{0})_{ν} v_{ν}^{μ} \frac{1}{\sqrt{n}})^{2} \frac{d_{0}^{2}}{2 η γ / a_{μ} + D + (n - 1) d}$	(3.7)
	$=$	$\frac{1}{p} p \sum μ = 1 (β^{0} \cdot v^{μ})^{2} \frac{d_{0}^{2} a_{μ}}{2 η γ / a_{μ} + D - d} + O (n)$
	$=$	$d_{0}^{2} ⟨ \frac{a^{2} (β^{0} \cdot v)^{2}}{2 η γ + (D - d) a} ⟩ + O (n)$

The averages in (3.6,3.7

) are now over the joint distribution of eigenvalues and eigenvectors of

A

only. Inserting the above RS expressions into (2.18), and using

C_{00} = {~ S}^{2}

, then gives us, with the short-hand

Dz=(2π)−12\rme−12z2\rmdz

$\frac{1}{n} Ψ_{R S} (\dots)$	$=$	$- \frac{1}{2} ζ (2 d_{0} c_{0} + D C - d c) + \frac{1}{2} (1 - ζ) log (2 π) - η ζ γ S^{2} + O (n)$	(3.8)
		$+ \frac{1}{2} [log (C - c) + \frac{c - c_{0}^{2} / C_{00}}{C - c}] - \frac{1}{2} ζ d_{0}^{2} ⟨ \frac{a^{2} (β^{0} \cdot v)^{2}}{2 η γ + (D - d) a} ⟩$
		$+ \frac{1}{2} ζ ⟨ log (\frac{2 η γ}{a} + D - d) ⟩ + \frac{1}{2} ζ ⟨ \frac{d a}{2 η γ + (D - d) a} ⟩ + \frac{1}{n} log ~ S$
		$−1nlog∫Dz∫\rmdy0√2π\rme−12B00y20∫\rmdt p(t\|y0,λ0)$
		$×[∫\rmdy \rme−12(B−b)y2+y(\rmi√b−b0y0)pγ(t\|y,λ)pγ(t\|y0,λ0)]n$
	$=$	$- \frac{1}{2} ζ (2 d_{0} c_{0} + D C - d c) + \frac{1}{2} (1 - ζ) log (2 π) - η ζ γ S^{2} + O (n)$
		$+ \frac{1}{2} [log (C - c) + \frac{c - c_{0}^{2} / {~ S}^{2}}{C - c}] - \frac{1}{2} ζ d_{0}^{2} ⟨ \frac{a^{2} (β^{0} \cdot v)^{2}}{2 η γ + (D - d) a} ⟩$
		$+ \frac{1}{2} ζ ⟨ log (\frac{2 η γ}{a} + D - d) ⟩ + \frac{1}{2} ζ ⟨ \frac{d a}{2 η γ + (D - d) a} ⟩ + \frac{1}{2 n} log ({~ S}^{2} B_{00})$
		$−1nlog∫DzDy0∫\rmdt p(t\|y0/√B00,λ0)$
		$×[∫\rmdy \rme−12(B−b)y2+y(\rmi√b−b0y0/√B00)pγ(t\|y,λ)pγ(t\|y0/√B00,λ0)]n$

We note that

	$B_{00}^{- 1} = {~ S}^{2} - n c_{0}^{2} / (C - c) + O (n^{2}), B - b = 1 / (C - c)$		(3.9)
	$b_{0} = - c_{0} / {~ S}^{2} (C - c) + O (n), b = \frac{c_{0}^{2} - c {~ S}^{2}}{{~ S}^{2} (C - c)^{2}} + O (n)$		(3.10)

and these identities enable us, after some simple rearrangements, to write the limit $Ψ_{R S} (\dots) = {lim}_{n \to 0} n^{- 1} Ψ_{R S} (\dots)$ in the much simpler form

$Ψ_{R S} (\dots)$	$=$	$- \frac{1}{2} ζ {2 d_{0} c_{0} + D C - d c + log (2 π) + 2 η γ S^{2}$	(3.11)
		$+ d_{0}^{2} ⟨ \frac{a^{2} (β^{0} \cdot v)^{2}}{2 η γ + (D - d) a} ⟩ - ⟨ log (\frac{2 η γ}{a} + D - d) ⟩ - ⟨ \frac{d a}{2 η γ + (D - d) a} ⟩}$
		$−∫DzDy0∫\rmdt p(t\|~Sy0,λ0)log∫Dy pγ$

Analysis of overfitting in the regularized Cox model

Authors

Replica analysis of overfitting in generalized linear models

Correction of overfitting bias in regression models

Structure Learning in Inverse Ising Problems Using ℓ_2-Regularized Linear Estimator

Optimal regularizations for data generation with probabilistic graphical models

Complex Structure Leads to Overfitting: A Structure Regularization Decoding Method for Natural Language Processing

A new analytical approach to consistency and overfitting in regularized empirical risk minimization

Towards the interpretation of time-varying regularization parameters in streaming penalized regression models

1 Introduction

2 Replica analysis of regularized Cox regression

2.1 Generalized replica formalism to include priors

2.2 Application to the regularized Cox proportional hazards model

2.3 Conversion into a saddle point problem

3 Replica symmetric theory

3.1 Replica symmetric saddle points

Analysis of overfitting in the regularized Cox model

Authors

Related Research

Replica analysis of overfitting in generalized linear models

Correction of overfitting bias in regression models

Structure Learning in Inverse Ising Problems Using ℓ_2-Regularized Linear Estimator

Optimal regularizations for data generation with probabilistic graphical models

Complex Structure Leads to Overfitting: A Structure Regularization Decoding Method for Natural Language Processing

A new analytical approach to consistency and overfitting in regularized empirical risk minimization

Towards the interpretation of time-varying regularization parameters in streaming penalized regression models

This week in AI

1 Introduction

2 Replica analysis of regularized Cox regression

2.1 Generalized replica formalism to include priors

2.2 Application to the regularized Cox proportional hazards model

2.3 Conversion into a saddle point problem

3 Replica symmetric theory

3.1 Replica symmetric saddle points