format: text

text:
# General Loss Bounds for Universal Sequence Prediction

## 1 Introduction

### 1.1 Induction

Many problems are of induction type, in which statements about the future have to be made, based on past observations. What is the probability of rain tomorrow, given the weather observations of the last few days? Is the Dow Jones likely to rise tomorrow, given the chart of the last years and possibly additional newspaper information? Can we reasonably doubt that the sun will rise tomorrow? Indeed, one definition of science is to predict the future, where, as an intermediate step, one tries to understand the past by developing theories and, as a consequence of prediction, one tries to manipulate the future. All induction problems may be studied in the Bayesian framework. The probability of observing $`x_t`$ at time $`t`$, given the observations $`x_1\mathrm{}x_{t1}`$ can be computed with Bayes’ rule, if we know the true probability distribution of observation sequences $`x_1x_2x_3\mathrm{}`$. The problem is that in many cases we do not even have a reasonable guess of the true distribution $`\mu `$. What is the true probability of weather sequences, stock charts, or sunrises?

### 1.2 Universal Sequence Prediction

Solomonoff \[Sol64\] had the idea to define a universal probability distribution111We use the term distribution slightly unprecisely for a probability measure. $`\xi `$ as a weighted average over all possible computable probability distributions. Lower weights were assigned to more complex distributions. He unified Epicurus’ principle of multiple explanations, Occams’ razor, and Bayes’ rule into an elegant formal theory. For a binary alphabet, the universal conditional probability used for predicting $`x_t`$ converges to the true conditional probability for $`t\mathrm{}`$ with probability 1. The convergence serves as a justification of using $`\xi `$ as a substitution for the usually unknown $`\mu `$. The framework can easily be generalized to other probability classes and weights \[Sol78\].

### 1.3 Contents

The main aim of this work is to prove expected loss bounds for general loss functions which measure the performance of $`\xi `$ relative to $`\mu `$, and to apply the results to games of chance. Details and proofs can be found in \[Hut01\]. There are good introductions and surveys of Solomonoff sequence prediction \[LV97\], inductive inference \[AS83, Sol97\], reasoning under uncertainty \[Grü98\], and competitive online statistics \[Vov99\] with interesting relations to this work. See \[Hut01\] and subsection 5.4 for details.

Section 2 explains notation and defines the generalized universal distribution $`\xi `$ as the $`w_{\mu _i}`$ weighted sum of probability distributions $`\mu _i`$ of a set $`M`$, which must include the true distribution $`\mu `$. This generalization is straightforward and causes no problems. $`\xi `$ multiplicatively dominates all $`\mu _iM`$, and the relative entropy between $`\mu `$ and $`\xi `$ is bounded by $`\mathrm{ln}\frac{1}{w_\mu }`$. Convergence of $`\xi `$ to $`\mu `$ is shown in Theorem 1.

Section 3 considers the case where a prediction or action $`y_t𝒴`$ results in a loss $`l_{x_ty_t}`$ if $`x_t`$ is the next symbol of the sequence. Optimal universal $`\mathrm{\Lambda }_\xi `$ and optimal informed $`\mathrm{\Lambda }_\mu `$ prediction schemes are defined for this case and loss bounds are proved. Theorems 2 and 3 bound the total loss $`L_\xi `$ of $`\mathrm{\Lambda }_\xi `$ by the total loss $`L_\mu `$ of $`\mathrm{\Lambda }_\mu `$ plus $`O(\sqrt{L_\mu })`$ terms.

Section 4 applies Theorem 3 to games of chance, defined as a sequence of bets, observations, and rewards. The average profit $`\overline{p}_{n\mathrm{\Lambda }_\xi }`$ achieved by the $`\mathrm{\Lambda }_\xi `$ scheme rapidly converges to the best possible average profit $`\overline{p}_{n\mathrm{\Lambda }_\mu }`$ achieved by the $`\mathrm{\Lambda }_\mu `$ scheme ($`\overline{p}_{n\mathrm{\Lambda }_\xi }\overline{p}_{n\mathrm{\Lambda }_\mu }=O(n^{1/2})`$). If there is a profitable scheme at all, asymptotically the universal $`\mathrm{\Lambda }_\xi `$ scheme will also become profitable. Theorem 4 lower bounds the time needed to reach the winning zone in terms of the relative entropy of $`\mu `$ and $`\xi `$. An attempt is made to give an information theoretic interpretation of the result.

Section 5 outlines possible extensions of the presented theory and results. They include arbitrary alphabets, partial, delayed and probabilistic prediction, classification, even more general loss functions, active systems influencing the environment, learning aspects, and a comparison to the weighted majority algorithm(s) and loss bounds.

## 2 Setup and Convergence

### 2.1 Strings and Probability Distributions

We denote binary strings by $`x_1x_2\mathrm{}x_n`$ with $`x_t\{0,1\}`$. We further use the abbreviations $`x_{n:m}:=x_nx_{n+1}\mathrm{}x_{m1}x_m`$ and $`x_{<n}:=x_1\mathrm{}x_{n1}`$. We use Greek letters for probability distributions. Let $`\rho (x_{1:t})`$ be the probability that an (infinite) sequence starts with $`x_1\mathrm{}x_t`$. The conditional probability

$$\rho (x_t|x_{<t})=\frac{\rho (x_{1:t})}{\rho (x_{<t})}$$
(1)
that a given string $`x_1\mathrm{}x_{t1}`$ is continued by $`x_t`$ is obtained by using Bayes’ rule. The prediction schemes will be based on these posteriors.

### 2.2 Universal Prior Probability Distribution

Every inductive inference problem can be brought into the following form: Given a string $`x_{<t}`$, take a guess at its continuation $`x_t`$. We will assume that the strings which have to be continued are drawn from a probability222This includes deterministic environments, in which case the probability distribution $`\mu `$ is $`1`$ for some sequence $`x_{1:\mathrm{}}`$ and $`0`$ for all others. We call probability distributions of this kind deterministic. distribution $`\mu `$. The maximal prior information a prediction algorithm can possess is the exact knowledge of $`\mu `$, but in many cases the true distribution is not known. Instead, the prediction is based on a guess $`\rho `$ of $`\mu `$. We expect that a predictor based on $`\rho `$ performs well, if $`\rho `$ is close to $`\mu `$ or converges, in a sense, to $`\mu `$. Let $`M:=\{\mu _1,\mu _2,\mathrm{}\}`$ be a finite or countable set of candidate probability distributions on strings. We define a weighted average on $`M`$

$$\begin{array}{c}\xi (x_{1:n}):=\underset{\mu _iM}{}w_{\mu _i}\mu _i(x_{1:n}),\\ \underset{\mu _iM}{}w_{\mu _i}=1,w_{\mu _i}>0.\end{array}$$
(2)
It is easy to see that $`\xi `$ is a probability distribution as the weights $`w_{\mu _i}`$ are positive and normalized to 1 and the $`\mu _iM`$ are probabilities. For finite $`M`$ a possible choice for the $`w`$ is to give all $`\mu _i`$ equal weight ($`w_{\mu _i}=\frac{1}{|M|}`$). We call $`\xi `$ universal relative to $`M`$, as it multiplicatively dominates all distributions in $`M`$

$$\xi (x_{1:n})w_{\mu _i}\mu _i(x_{1:n})\text{for all}\mu _iM.$$
(3)
In the following, we assume that $`M`$ is known and contains the true distribution, i.e. $`\mu M`$. This is not a serious constraint if we include all computable probability distributions in $`M`$ with a high weight assigned to simple $`\mu _i`$. Solomonoff’s universal semi-measure is obtained if we include all enumerable semi-measures in $`M`$ with weights $`w_{\mu _i}2^{K(\mu _i)}`$, where $`K(\mu _i)`$ is the length of the shortest program for $`\mu _i`$ \[Sol64, Sol78, LV97\]. A detailed discussion of various general purpose choices for $`M`$ is given in \[Hut01\].

Furthermore, we need the relative entropy between $`\mu `$ and $`\xi `$:

$$h_t(x_{<t}):=\underset{x_t\{0,1\}}{}\mu (x_t|x_{<t})\mathrm{ln}\frac{\mu (x_t|x_{<t})}{\xi (x_t|x_{<t})}$$
(4)
$`H_n`$ is then defined as the sum-expectation, for which the following can be shown

$$H_n:=\underset{t=1}{\overset{n}{}}\underset{x_{<t}\{0,1\}^{t1}}{}\mu (x_{<t})h_t(x_{<t})\mathrm{ln}\frac{1}{w_\mu }=:d_\mu $$
(5)
The following theorem shows the important property of $`\xi `$ converging to the true distribution $`\mu `$, in a sense.

###### Theorem 1 (Convergence)

Let there be binary sequences $`x_1x_2\mathrm{}`$ drawn with probability $`\mu (x_{1:n})`$ for the first $`n`$ symbols. The universal conditional probability $`\xi (x_t|x_{<t})`$ of the next symbol $`x_t`$ given $`x_{<t}`$ is related to the true conditional probability $`\mu (x_t|x_{<t})`$ in the following way:

$$\begin{array}{cc}\hfill i)& \underset{t=1}{\overset{n}{}}\underset{x_{1:t}\{0,1\}^t}{}\mu (x_{<t})\left(\mu (x_t|x_{<t})\xi (x_t|x_{<t})\right)^2\hfill \\ & H_nd_\mu =\mathrm{ln}\frac{1}{w_\mu }<\mathrm{}\\ \hfill ii)& \xi (x_t|x_{<t})\mu (x_t|x_{<t})\text{for }t\mathrm{}\text{ with}\hfill \\ & \mu \text{ probability 1}\end{array}$$

where $`H_n`$ is the relative entropy (5), and $`w_\mu `$ is the weight (2) of $`\mu `$ in $`\xi `$.

$`(i)`$ and (5) are easy generalizations of \[Sol78\] to arbitrary weights $`w_\mu `$ and an arbitrary probability set $`M`$. For $`n\mathrm{}`$ the l.h.s. of $`(i)`$ is an infinite $`t`$-sum over positive arguments, which is bounded by the finite constant $`d_\mu `$ on the r.h.s. Hence the arguments must converge to zero for $`t\mathrm{}`$. Since the arguments are $`\mu `$ expectations of the squared difference of $`\xi `$ and $`\mu `$, this means that $`\xi (x_t|x_{<t})`$ converges333More precisely $`\xi (x_t|x_{<t})\mu (x_t|x_{<t})`$ converges to zero for $`t\mathrm{}`$ with $`\mu `$ probability 1 or, more stringent, in a mean squared sense. to $`\mu (x_t|x_{<t})`$ with $`\mu `$ probability 1. This proves $`(ii)`$. Since the conditional probabilities are the basis of all prediction algorithms considered in this work, we expect a good prediction performance if we use $`\xi `$ as a guess of $`\mu `$. Performance measures are defined in the next section.

## 3 Loss Bounds

### 3.1 Unit Loss Function

A prediction is very often the basis for some decision. The decision results in an action, which itself leads to some reward or loss. If the action itself can influence the environment we enter the domain of acting agents which has been analyzed in the context of universal probability in \[Hut00\]. To stay in the framework of (passive) prediction we have to assume that the action itself does not influence the environment. Let $`l_{x_ty_t}IR`$ be the received loss when taking action $`y_t𝒴`$ and $`x_t\{0,1\}`$ is the $`t^{th}`$ symbol of the sequence. We demand $`l`$ to be normalized, i.e. $`0l_{x_ty_t}1`$. For instance, if we make a sequence of weather forecasts $`\{0,1\}=\{`$sunny, rainy$`\}`$ and base our decision, whether to take an umbrella or wear sunglasses $`𝒴=\{`$umbrella, sunglasses$`\}`$ on it, the action of taking the umbrella or wearing sunglasses does not influence the future weather (ignoring the butterfly effect). Reasonable losses may be

| Loss | sunny | rainy |
| --- | --- | --- |
| umbrella | 0.3 | 0.1 |
| sunglasses | 0.0 | 1.0 |

In many cases the prediction of $`x_t`$ can be identified or is already the action $`y_t`$. The forecast sunny can be identified with the action wear sunglasses, and rainy with take umbrella. In the following, we assume “predictive” actions of this kind, i.e. $`𝒴=\{0,1\}`$. General action spaces $`𝒴`$ and general alphabets $`𝒜`$ are considered in \[Hut01\].

The true probability of the next symbol being $`x_t`$, given $`x_{<t}`$, is $`\mu (x_t|x_{<t})`$. The expected loss when predicting $`y_t`$ is $`\mu (1|x_{<t})l_{1y_t}+\mu (0|x_{<t})l_{0y_t}`$. The goal is to minimize the expected loss. More generally we define the $`\mathrm{\Lambda }_\rho `$ prediction scheme

$$y_t^{\mathrm{\Lambda }_\rho }:=\text{arg}\underset{y_t}{\mathrm{min}}\underset{x_t\{0,1\}}{}\rho (x_t|x_{<t})l_{x_ty_t}$$
(6)
which minimizes the $`\rho `$-expected loss. This is a threshold strategy with $`y_t^{\mathrm{\Lambda }_\rho }=0/1`$ for $`\rho (1|x_{<t})_<^>\gamma `$, where $`\gamma :=\frac{l_{01}l_{00}}{l_{01}l_{00}+l_{10}l_{11}}`$. As the true distribution is $`\mu `$, the actual $`\mu `$ expected loss when $`\mathrm{\Lambda }_\rho `$ predicts the $`t^{th}`$ symbol and the total $`\mu `$-expected loss in the first $`n`$ predictions are

$$\begin{array}{c}l_{t\mathrm{\Lambda }_\rho }(x_{<t}):=\underset{x_t}{}\mu (x_t|x_{<t})l_{x_ty_t^{\mathrm{\Lambda }_\rho }},\\ L_{n\mathrm{\Lambda }_\rho }:=\underset{t=1}{\overset{n}{}}\underset{x_{<t}}{}\mu (x_{<t})l_{t\mathrm{\Lambda }_\rho }(x_{<t}).\end{array}$$
(7)
In the special case $`l_{01}=l_{10}=1`$ and $`l_{00}=l_{11}=0`$, the bit with the highest $`\rho `$ probability is predicted ($`\gamma =\frac{1}{2}`$), and $`L_{n\mathrm{\Lambda }_\rho }`$ is the total expected number of prediction errors.

If $`\mu `$ is known, $`\mathrm{\Lambda }_\mu `$ is obviously the best prediction scheme in the sense of achieving minimal expected loss

$$L_{n\mathrm{\Lambda }_\mu }L_{n\mathrm{\Lambda }_\rho }\text{for any}\mathrm{\Lambda }_\rho $$
(8)
The predictor $`\mathrm{\Lambda }_\xi `$, based on the universal distribution $`\xi `$, is of special interest.

###### Theorem 2 (Unit loss bound)

Let there be binary sequences $`x_1x_2\mathrm{}`$ drawn with probability $`\mu (x_{1:n})`$ for the first $`n`$ symbols. A system predicting $`y_t\{0,1\}`$ given $`x_{<t}`$ receives loss $`l_{x_ty_t}[0,1]`$ if $`x_t`$ is the true $`t^{th}`$ symbol of the sequence. The $`\mathrm{\Lambda }_\rho `$-system (6) predicts as to minimize the $`\rho `$-expected loss. $`\mathrm{\Lambda }_\xi `$ is the universal prediction scheme based on the universal prior $`\xi `$. $`\mathrm{\Lambda }_\mu `$ is the optimal informed prediction scheme. The total $`\mu `$-expected losses $`L_{n\mathrm{\Lambda }_\xi }`$ of $`\mathrm{\Lambda }_\xi `$ and $`L_{n\mathrm{\Lambda }_\mu }`$ of $`\mathrm{\Lambda }_\mu `$ as defined in (7) are bounded in the following way

$$\begin{array}{c}0L_{n\mathrm{\Lambda }_\xi }L_{n\mathrm{\Lambda }_\mu }H_n+\sqrt{4L_{n\mathrm{\Lambda }_\mu }H_n+H_n^2}\end{array}$$

where $`H_n\mathrm{ln}\frac{1}{w_\mu }`$ is the relative entropy (5), and $`w_\mu `$ is the weight (2) of $`\mu `$ in $`\xi `$.

First, we observe that the total loss $`L_{\mathrm{}\mathrm{\Lambda }_\xi }`$ of the universal $`\mathrm{\Lambda }_\xi `$ predictor is finite if the total loss $`L_{\mathrm{}\mathrm{\Lambda }_\mu }`$ of the informed $`\mathrm{\Lambda }_\mu `$ predictor is finite. This is especially the case for deterministic $`\mu `$ and $`l_{00}=l_{11}=0`$, as $`L_{n\mathrm{\Lambda }_\mu }0`$ in this case444Remember that we named a probability distribution deterministic if it is 1 for exactly one sequence and 0 for all others., i.e. $`\mathrm{\Lambda }_\xi `$ receives a finite loss on deterministic environments if a correct prediction results in zero loss. More precisely, $`L_{\mathrm{}\mathrm{\Lambda }_\xi }2H_{\mathrm{}}2\mathrm{ln}\frac{1}{w_\mu }`$. A combinatoric argument shows that there are $`M`$ and $`\mu M`$ with $`L_{\mathrm{}\mathrm{\Lambda }_\xi }\mathrm{log}_2|M|`$. This shows that the upper bound $`L_{\mathrm{}\mathrm{\Lambda }_\xi }2\mathrm{ln}|M|`$ for uniform $`w`$ is rather tight. For more complicated probabilistic environments, where even the ideal informed system makes an infinite number of errors, the theorem ensures that the loss excess $`L_{n\mathrm{\Lambda }_\xi }L_{n\mathrm{\Lambda }_\mu }`$ is only of order $`\sqrt{L_{n\mathrm{\Lambda }_\mu }}`$. The excess is quantified in terms of the information content $`H_n`$ of $`\mu `$ (relative to $`\xi `$), or the weight $`w_\mu `$ of $`\mu `$ in $`\xi `$. This ensures that the loss densities $`L_n/n`$ of both systems converge to each other for $`n\mathrm{}`$. Actually, the theorem ensures more, namely that the quotient converges to 1, and also gives the speed of convergence $`L_{n\mathrm{\Lambda }_\xi }/L_{n\mathrm{\Lambda }_\mu }=1+O(L_{n\mathrm{\Lambda }_\mu }^{1/2})1`$ for $`L_{n\mathrm{\Lambda }_\mu }\mathrm{}`$.

### 3.2 Proof Sketch of Theorem 2

The first inequality in Theorem 2 has already been proved (8). For the second inequality, let us start more modestly and try to find constants $`A>0`$ and $`B>0`$ that satisfy the linear inequality

$$L_{n\mathrm{\Lambda }_\xi }(A+1)L_{n\mathrm{\Lambda }_\mu }+(B+1)H_n.$$
(9)
If we could show

$$l_{t\mathrm{\Lambda }_\xi }(x_{<t})A^{}l_{t\mathrm{\Lambda }_\mu }(x_{<t})+B^{}h_t(x_{<t})$$
(10)
with $`A^{}:=A+1`$ and $`B^{}:=B+1`$ for all $`tn`$ and all $`x_{<t}`$, (9) would follow immediately by summation and the definition of $`L_n`$ and $`H_n`$. With the abbreviations

$$i=x_t,y_i=\mu (x_t|x_{<t}),z_i=\xi (x_t|x_{<t})$$

$$m=y_t^{\mathrm{\Lambda }_\mu },s=y_t^{\mathrm{\Lambda }_\xi }$$

the loss and entropy can be expressed by $`l_{t\mathrm{\Lambda }_\xi }=_iy_il_{is}`$, $`l_{t\mathrm{\Lambda }_\mu }=_iy_il_{im}`$ and $`h_t=_iy_i\mathrm{ln}\frac{y_i}{z_i}`$. Inserting this into (10) and rearranging terms we have to prove

$$B^{}\underset{i=0}{\overset{1}{}}y_i\mathrm{ln}\frac{y_i}{z_i}+\underset{i=0}{\overset{1}{}}y_i(A^{}l_{im}l_{is})\stackrel{\mathrm{?}}{}\mathrm{\hspace{0.33em}0}.$$
(11)
By definition (6) of $`y_t^{\mathrm{\Lambda }_\mu }`$ and $`y_t^{\mathrm{\Lambda }_\xi }`$ we have

$$\underset{i}{}y_il_{im}\underset{i}{}y_il_{ij}\text{and}\underset{i}{}z_il_{is}\underset{i}{}z_il_{ij}$$
(12)
for all $`j`$. Actually, we need the first constraint only for $`j=s`$ and the second for $`j=m`$. The cases $`l_{im}>l_{is}i`$ and $`l_{is}>l_{im}i`$ contradict the first/second inequality (12). Hence we can assume $`l_{0m}l_{0s}`$ and $`l_{1m}l_{1s}`$. The symmetric case $`l_{0m}l_{0s}`$ and $`l_{1m}l_{1s}`$ is proved analogously or can be reduced to the first case by renumbering the indices ($`01`$). Using the abbreviations $`a:=l_{0m}l_{0s}`$, $`b:=l_{1s}l_{1m}`$, $`c:=y_1l_{1m}+y_0l_{0s}`$, $`y=y_1=1y_0`$ and $`z=z_1=1z_0`$ we can write (11) as

$$f(y,z):=$$
(13)
$$B^{}[y\mathrm{ln}\frac{y}{z}+(1y)\mathrm{ln}\frac{1y}{1z}]+A^{}(1y)ayb+Ac\stackrel{\mathrm{?}}{}\mathrm{\hspace{0.33em}0}$$

for $`zb(1z)a`$ and $`0a,b,c,y,z1`$. The constraint (12) on $`y`$ has been dropped since (13) will turn out to be true for all $`y`$. Furthermore, we can assume that $`d:=A^{}(1y)ayb0`$ since for $`d>0`$, $`f`$ is trivially positive ($`h_t0`$). Multiplying $`d`$ with a constant $`1`$ will decrease $`f`$. Let us first consider the case $`z\frac{1}{2}`$. We multiply the $`d`$ term by $`1/b1`$, i.e. replace it with $`A^{}(1y)\frac{a}{b}y`$. From the constraint on $`z`$ we known that $`\frac{a}{b}\frac{z}{1z}`$. We can decrease $`f`$ further by replacing $`\frac{a}{b}`$ by $`\frac{z}{1z}`$ and by dropping $`Ac`$. Hence, (13) is proved for $`z\frac{1}{2}`$ if we can prove

$$B^{}[\mathrm{}]+A^{}(1y)\frac{z}{1z}y\stackrel{\mathrm{?}}{}\mathrm{\hspace{0.33em}0}\text{for}z\frac{1}{2}.$$
(14)
The case $`z\frac{1}{2}`$ is treated similarly. We scale $`d`$ with $`1/a1`$, i.e. replace it with $`A^{}(1y)y\frac{b}{a}`$. From the constraint on $`z`$ we know that $`\frac{b}{a}\frac{1z}{z}`$. We decrease $`f`$ further by replacing $`\frac{b}{a}`$ by $`\frac{1z}{z}`$ and by dropping $`Ac`$. Hence (13) is proved for $`z\frac{1}{2}`$ if we can prove

$$B^{}[\mathrm{}]+A^{}(1y)y\frac{1z}{z}\stackrel{\mathrm{?}}{}\mathrm{\hspace{0.33em}0}\text{for}z\frac{1}{2}.$$
(15)
In \[Hut01\] we prove that (14) and (15) indeed hold for $`B\frac{1}{4}A+\frac{1}{A}`$. The cautious reader may check the inequalities numerically. So in summary we proved that (9) holds for $`B\frac{1}{4}A+\frac{1}{A}`$. Inserting $`B=\frac{1}{4}A+\frac{1}{A}`$ into (9) and minimizing the r.h.s. with respect to $`A`$ leads to the bound of Theorem 2 (with $`A^2=H_n/(L_{n\mathrm{\Lambda }_\mu }+\frac{1}{4}H_n)`$) $``$.

### 3.3 General Loss

There are only very few restrictions imposed on the loss $`l_{x_ty_t}`$ in Theorem 2, namely that it is static and in the unit interval $`[0,1]`$. If we look at the proof of Theorem 2, we see that the time-independence has not been used at all. The proof is still valid for an individual loss function $`l_{x_ty_t}^t[0,1]`$ for each step $`t`$. The loss might even depend on the actual history $`x_{<t}`$. The case of a loss $`l_{x_ty_t}^t(x_{<t})`$ bounded to a general interval $`[l_{min},l_{max}]`$ can be reduced to the unit interval case by rescaling $`l`$. We introduce a scaled loss $`l^{}`$

$$0l_{}^{}{}_{x_ty_t}{}^{t}(x_{<t}):=\frac{l_{x_ty_t}^t(x_{<t})l_{min}}{l_\mathrm{\Delta }}\mathrm{\hspace{0.33em}1},$$

$$\text{where}l_\mathrm{\Delta }:=l_{max}l_{min}.$$

The prediction scheme $`\mathrm{\Lambda }_\rho ^{}`$ based on $`l^{}`$ is identical to the original prediction scheme $`\mathrm{\Lambda }_\rho `$ based on $`l`$, since arg$`\mathrm{min}`$ in (6) is not affected by a constant scaling and a shift of its argument. From $`y_t^{\mathrm{\Lambda }_\rho ^{}}=y_t^{\mathrm{\Lambda }_\rho }`$ it follows that $`l_{t\mathrm{\Lambda }_\rho }^{}=(l_{t\mathrm{\Lambda }_\rho }l_{min})/l_\mathrm{\Delta }`$ and $`L_{n\mathrm{\Lambda }_\rho }^{}=(L_{n\mathrm{\Lambda }_\rho }l_{min})/l_\mathrm{\Delta }`$ ($`H_n^{}H_n`$, since $`l`$ is not involved). Theorem 2 is valid for the primed quantities, since $`l^{}[0,1]`$. Inserting $`L_{n\mathrm{\Lambda }_{\mu /\xi }}^{}`$ and rearranging terms we get

###### Theorem 3 (General loss bound)

Let there be binary sequences $`x_1x_2\mathrm{}`$ drawn with probability $`\mu (x_{1:n})`$ for the first $`n`$ symbols. A system taking action (or predicting) $`y_t𝒴`$ given $`x_{<t}`$ receives loss $`l_{x_ty_t}^t(x_{<t})[l_{min},l_{min}+l_\mathrm{\Delta }]`$ if $`x_t`$ is the true $`t^{th}`$ symbol of the sequence. The $`\mathrm{\Lambda }_\rho `$-system (6) acts (or predicts) as to minimize the $`\rho `$-expected loss. $`\mathrm{\Lambda }_\xi `$ is the universal prediction scheme based on the universal prior $`\xi `$. $`\mathrm{\Lambda }_\mu `$ is the optimal informed prediction scheme. The total $`\mu `$-expected losses $`L_{n\mathrm{\Lambda }_\xi }`$ and $`L_{n\mathrm{\Lambda }_\mu }`$ of $`\mathrm{\Lambda }_\xi `$ and $`\mathrm{\Lambda }_\mu `$ as defined in (7) are bounded in the following way

$$0L_{n\mathrm{\Lambda }_\xi }L_{n\mathrm{\Lambda }_\mu }$$

$$l_\mathrm{\Delta }H_n+\sqrt{4(L_{n\mathrm{\Lambda }_\mu }nl_{min})l_\mathrm{\Delta }H_n+l_\mathrm{\Delta }^2H_n^2}$$

where $`H_n\mathrm{ln}\frac{1}{w_\mu }`$ is the relative entropy (5), and $`w_\mu `$ is the weight (2) of $`\mu `$ in $`\xi `$.

## 4 Application to Games of Chance

### 4.1 Introduction/Example

Think of investing in the stock market. At time $`t`$ an amount of money $`s_t`$ is invested in portfolio $`y_t`$, where we have access to past knowledge $`x_{<t}`$ (e.g. charts). After our choice of investment we receive new information $`x_t`$, and the new portfolio value is $`r_t`$. The best we can expect is to have a probabilistic model $`\mu `$ of the behaviour of the stock-market. The goal is to maximize the net $`\mu `$-expected profit $`p_t=r_ts_t`$. Nobody knows $`\mu `$, but the assumption of all traders is that there is a computable, profitable $`\mu `$ they try to find or approximate. From Theorem 1 we know that Solomonoff’s universal prior $`\xi (x_t|x_{<t})`$ converges to any computable $`\mu (x_t|x_{<t})`$ with probability 1. If there is a computable, asymptotically profitable trading scheme at all, the $`\mathrm{\Lambda }_\xi `$ scheme should also be profitable in the long run. To get a practically useful, computable scheme we have to restrict $`M`$ to a finite set of computable distributions, e.g. with bounded Levin complexity $`Kt`$ \[LV97\]. Although convergence of $`\xi `$ to $`\mu `$ is pleasing, what we are really interested in is whether $`\mathrm{\Lambda }_\xi `$ is asymptotically profitable and how long it takes to become profitable. This will be explored in the following.

### 4.2 Games of Chance

We use Theorem 3 (or its generalization to arbitrary action and alphabet, proved in \[Hut01\]) to estimate the time needed to reach the winning threshold when using $`\mathrm{\Lambda }_\xi `$ in a game of chance. We assume a game (or a sequence of possibly correlated games) which allows a sequence of bets and observations. In step $`t`$ we bet, depending on the history $`x_{<t}`$, a certain amount of money $`s_t`$, take some action $`y_t`$, observe outcome $`x_t`$, and receive reward $`r_t`$. Our profit, which we want to maximize, is $`p_t=r_ts_t`$. The loss, which we want to minimize, can be defined as the negative profit, $`l_{x_ty_t}=p_t`$. The probability of outcome $`x_t`$, possibly depending on the history $`x_{<t}`$, is $`\mu (x_t|x_{<t})`$. The total $`\mu `$ expected profit when using scheme $`\mathrm{\Lambda }_\rho `$ is $`P_{n\mathrm{\Lambda }_\rho }=L_{n\mathrm{\Lambda }_\rho }`$. If we knew $`\mu `$, the optimal strategy to maximize our expected profit is just $`\mathrm{\Lambda }_\mu `$. We assume $`P_{n\mathrm{\Lambda }_\mu }>0`$ (otherwise there is no winning strategy at all, since $`P_{n\mathrm{\Lambda }_\mu }P_{n\mathrm{\Lambda }_\rho }\rho `$). Often we are not in the favorable position of knowing $`\mu `$, but we know (or assume) that $`\mu M`$ for some $`M`$, for instance that $`\mu `$ is a computable probability distribution. From Theorem 3 we see that the average profit per round $`\overline{p}_{n\mathrm{\Lambda }_\xi }:=\frac{1}{n}P_{n\mathrm{\Lambda }_\xi }`$ of the universal $`\mathrm{\Lambda }_\xi `$ scheme converges to the average profit per round $`\overline{p}_{n\mathrm{\Lambda }_\mu }:=\frac{1}{n}P_{n\mathrm{\Lambda }_\mu }`$ of the optimal informed scheme, i.e. asymptotically we can make the same money even without knowing $`\mu `$, by just using the universal $`\mathrm{\Lambda }_\xi `$ scheme. Theorem 3 allows us to lower bound the universal profit $`P_{n\mathrm{\Lambda }_\xi }`$

$$P_{n\mathrm{\Lambda }_\xi }P_{n\mathrm{\Lambda }_\mu }p_\mathrm{\Delta }H_n\sqrt{4(np_{max}P_{n\mathrm{\Lambda }_\mu })p_\mathrm{\Delta }H_n+p_\mathrm{\Delta }^2H_n^2}$$
(16)
where $`p_{max}`$ is the maximal profit per round and $`p_\mathrm{\Delta }`$ the profit range. The time needed for $`\mathrm{\Lambda }_\xi `$ to perform well can also be estimated. An interesting quantity is the expected number of rounds needed to reach the winning zone. Using $`P_{n\mathrm{\Lambda }_\mu }>0`$ one can show that the r.h.s. of (16) is positive if, and only if

$$n>\frac{2p_\mathrm{\Delta }(2p_{max}\overline{p}_{n\mathrm{\Lambda }_\mu })}{\overline{p}_{n\mathrm{\Lambda }_\mu }^2}H_n.$$
(17)
###### Theorem 4 (Time to Win)

Let there be binary sequences $`x_1x_2\mathrm{}`$ drawn with probability $`\mu (x_{1:n})`$ for the first $`n`$ symbols. In step $`t`$ we make a bet, depending on the history $`x_{<t}`$, take some action $`y_t`$, and observe outcome $`x_t`$. Our net profit is $`p_t[p_{max}p_\mathrm{\Delta },p_{max}]`$. The $`\mathrm{\Lambda }_\rho `$-system (6) acts as to maximize the $`\rho `$-expected profit. $`P_{n\mathrm{\Lambda }_\rho }`$ is the total and $`\overline{p}_{n\mathrm{\Lambda }_\rho }=\frac{1}{n}P_{n\mathrm{\Lambda }_\rho }`$ is the average expected profit of the first $`n`$ rounds. For the universal $`\mathrm{\Lambda }_\xi `$ and for the optimal informed $`\mathrm{\Lambda }_\mu `$ prediction scheme the following holds:

$$\begin{array}{cc}\hfill i)& \overline{p}_{n\mathrm{\Lambda }_\xi }=\overline{p}_{n\mathrm{\Lambda }_\mu }O(n^{1/2})\overline{p}_{n\mathrm{\Lambda }_\mu }\text{for}n\mathrm{}\hfill \\ \hfill ii)& \text{if}n>\left(\frac{2p_\mathrm{\Delta }}{\overline{p}_{n\mathrm{\Lambda }_\mu }}\right)^2d_\mu \text{and}\overline{p}_{n\mathrm{\Lambda }_\mu }>0\overline{p}_{n\mathrm{\Lambda }_\xi }>0\hfill \end{array}$$

where $`w_\mu =e^{d_\mu }`$ is the weight (2) of $`\mu `$ in $`\xi `$.

By dividing (16) by $`n`$ and using $`H_nd_\mu `$ (5) we see that the leading order of $`\overline{p}_{n\mathrm{\Lambda }_\xi }\overline{p}_{n\mathrm{\Lambda }_\mu }`$ is bounded by $`\sqrt{4p_\mathrm{\Delta }p_{max}d_\mu /n}`$, which proves $`(i)`$. The condition in $`(ii)`$ is actually a weakening of (17). $`P_{n\mathrm{\Lambda }_\xi }`$ is trivially positive for $`p_{min}>0`$, since in this wonderful case all profits are positive. For negative $`p_{min}`$ the condition of $`(ii)`$ implies (17), since $`p_\mathrm{\Delta }>p_{max}`$, and (17) implies positive (16), i.e. $`P_{n\mathrm{\Lambda }_\xi }>0`$, which proves $`(ii)`$.

If a winning strategy $`\mathrm{\Lambda }_\rho `$ with $`\overline{p}_{n\mathrm{\Lambda }_\rho }>\epsilon >0`$ exists, then $`\mathrm{\Lambda }_\xi `$ is asymptotically also a winning strategy with the same average profit.

### 4.3 Information-Theoretic Interpretation

We try to give an intuitive explanation of Theorem 4$`(ii)`$. We know that $`\xi (x_t|x_{<t})`$ converges to $`\mu (x_t|x_{<t})`$ for $`t\mathrm{}`$. In a sense $`\mathrm{\Lambda }_\xi `$ learns $`\mu `$ from past data $`x_{<t}`$. The information content in $`\mu `$ relative to $`\xi `$ is $`\mathrm{ln}2H_{\mathrm{}}d_\mu \mathrm{ln}2`$. One might think of a Shannon-Fano prefix code of $`\mu _iM`$ of length $`{}_{}{}^{}d_{\mu _i}^{}\mathrm{ln}2^{}`$, which exists since the Kraft inequality $`_i2^{^{}d_{\mu _i}\mathrm{ln}2^{}}_iw_{\mu _i}1`$ is satisfied. $`d_\mu \mathrm{ln}2`$ bits have to be learned before $`\mathrm{\Lambda }_\xi `$ can be as good as $`\mathrm{\Lambda }_\mu `$. In the worst case, the only information contained in $`x_t`$ is in form of the received profit $`p_t`$. Remember that we always know the profit $`p_t`$ before the next cycle starts.

Assume that the distribution of the profits in the interval $`[p_{min},p_{max}]`$ is mainly due to noise, and there is only a small informative signal of amplitude $`\overline{p}_{n\mathrm{\Lambda }_\mu }`$. To reliably determine the sign of a signal of amplitude $`\overline{p}_{n\mathrm{\Lambda }_\mu }`$, disturbed by noise of amplitude $`p_\mathrm{\Delta }`$, we have to resubmit a bit $`O((p_\mathrm{\Delta }/\overline{p}_{n\mathrm{\Lambda }_\mu })^2)`$ times (this reduces the standard deviation below the signal amplitude $`\overline{p}_{n\mathrm{\Lambda }_\mu }`$). To learn $`\mu `$, $`d_\mu \mathrm{ln}2`$ bits have to be transmitted, which requires $`nO((p_\mathrm{\Delta }/\overline{p}_{n\mathrm{\Lambda }_\mu })^2)d_\mu \mathrm{ln}2`$ cycles. This expression coincides with the condition in $`(ii)`$. Identifying the signal amplitude with $`\overline{p}_{n\mathrm{\Lambda }_\mu }`$ is the weakest part of this consideration, as we have no argument why this should be true. It may be interesting to make the analogy more rigorous, which may also lead to a simpler proof of $`(ii)`$ not based on Theorems 2 and 3.

## 5 Outlook

In the following we discuss several directions in which the findings of this work may be extended.

### 5.1 General Alphabet

In many, cases the prediction unit is not a bit, but a letter from a finite alphabet $`𝒜`$. Non-binary prediction cannot be (easily) reduced to the binary case. One might think of a binary coding of the symbols $`x_t𝒜`$ in the sequence $`x_1x_2\mathrm{}`$. But this makes it necessary to predict a block of bits $`x_t`$, before one receives the true block of bits $`x_t`$, which differs from the bit by bit prediction, considered here and in \[Sol78\]! Fortunately, all theorems (1-4) take over to general alphabet \[Hut01\]. Unfortunately, the proofs are rather complex. In many cases the basic prediction unit is not even a letter from a finite alphabet, but a number (for inducing number sequences), or a word (for completing sentences), a real number or vector (for physical measurements). The prediction may either be generalized to a block by block prediction of symbols or, more suitably, the finite alphabet $`𝒜`$ could be generalized to countable (numbers, words) or continuous (real or vector) alphabet. The theorems should generalize to countably infinite alphabets by appropriately taking the limit $`|𝒜|\mathrm{}`$ and to continuous alphabets by a denseness or separability argument.

### 5.2 Partial Prediction, Delayed Prediction, Classification

The $`\mathrm{\Lambda }_\rho `$ schemes may also be used for partial prediction where, for instance, only every $`m^{th}`$ symbol is predicted. This can be arranged by setting the loss $`l^t`$ to zero when no prediction is made, e.g. if $`t`$ is not a multiple of $`m`$. Classification could be interpreted as partial sequence prediction, where $`x_{(t1)m+1:km1}`$ is classified as $`x_{km}`$. There are better ways for classification by treating $`x_{(t1)m+1:km1}`$ as pure conditions in $`\xi `$, as has been done in \[Hut00\] in a more general context. Another possibility is to generalize the prediction schemes and theorems to delayed sequence prediction, where the true symbol $`x_t`$ is given only in cycle $`t+d`$. A delayed feedback is common in many practical problems.

### 5.3 More Active Systems

Prediction means guessing the future, but not influencing it. We mentioned the possibility of interpreting $`y_t𝒴`$ as an action with $`𝒴𝒜`$. This tiny step towards a more active system is described in more detail in \[Hut01\]. The probability $`\mu `$ is still independent of the action, and the loss function $`l^t`$ has to be known in advance. This ensures that the greedy strategy (6) is optimal. The loss function may be generalized to depend not only on the history $`x_{<t}`$, but also on the historic actions $`y_{<t}`$ with $`\mu `$ still independent of the action. It would be interesting to know whether the scheme $`\mathrm{\Lambda }`$ and/or the loss bounds generalize to this case. The full model of an acting agent influencing the environment has been developed in \[Hut00\], but loss bounds have yet to be proven.

### 5.4 The Weighted Majority Algorithm(s)

The Weighted Majority (WM) algorithm is a related universal forecasting algorithm. It was invented by Littlestone and Warmuth \[LW89, LW94\] and Vovk \[Vov92\] and further developed in \[Ces97, HKW98, KW99\] and others. Many variations known by many names have meanwhile been invented. Early works in this direction are \[Daw84, Ris89\]. See \[Vov99\] for a review and further references. The setting and basic idea of WM are the following. Consider a finite binary sequence $`x_1x_2\mathrm{}x_n\{0,1\}^n`$ and a finite set $``$ of experts $`e`$ making predictions $`x_t^e`$ in the unit interval $`[0,1]`$ based on past observations $`x_1x_2\mathrm{}x_{t1}`$. The loss of expert $`e`$ in step $`t`$ is defined as $`|x_tx_t^e|`$. In the case of binary predictions $`x_t^e\{0,1\}`$, $`|x_tx_t^e|`$ coincides with our error measure defined in \[Hut01\]. The WM algorithm $`p_{\beta n}`$ combines the predictions of all experts. It forms its own prediction $`x_t^p`$ according to some weighted average of the expert’s predictions $`x_t^e`$. There are certain update rules for the weights depending on some parameter $`\beta `$. Various bounds for the total loss $`L_p(𝐱):=_{t=1}^n|x_tx_t^p|`$ of WM in terms of the total loss $`L_\epsilon (𝐱):=_{t=1}^n|x_tx_t^\epsilon |`$ of the best expert $`\epsilon `$ have been proven. It is possible to fine tune $`\beta `$ and to eliminate the necessity of knowing $`n`$ in advance. The most general bound of this kind is \[Ces97\]

$$L_p(𝐱)L_\epsilon (𝐱)+2.8\mathrm{ln}||+4\sqrt{L_\epsilon (𝐱)\mathrm{ln}||}.$$
(18)
It is interesting that our bound in Theorem 2 (with $`H_n\mathrm{ln}|M|`$ for uniform weights) has a quite similar structure as this bound, although the algorithms, the settings, the proofs and the interpretation are quite different. Whereas WM performs well in any environment, but only relative to a given set of experts $``$, our $`\mathrm{\Lambda }_\xi `$ predictor competes with the best possible $`\mathrm{\Lambda }_\mu `$ predictor (and hence with any other $`\rho `$ predictor), but only for a given set of environments $`M`$. WM depends on the set of expert, $`\mathrm{\Lambda }_\xi `$ depends on the set of environments $`M`$. The basic $`p_{\beta n}`$ algorithm has been extended in different directions: incorporation of different initial weights ($`||\mathrm{ln}\frac{1}{w_i}`$) \[LW89, Vov92\], more general loss functions \[HKW98\], continuous valued outcomes \[HKW98\], and multi-dimensional predictions \[KW99\] (but not yet for the absolute loss). The works of Yamanishi \[Yam97\] and \[Yam98\] lie somewhat in between WM and this work; “WM” techniques are used to prove expected loss bounds (but only for sequences of independent symbols/experiments and different classes of loss functions). Finally, note that the predictions of WM are continuous. In a sense it is more natural to predict $`0`$ or $`1`$ on a binary sequence, rather than some real number. On the other hand it is possible to convert the continuous prediction of WM into a probabilistic binary prediction by interpreting $`x_t^p[0,1]`$ as the probability of predicting $`1`$, and $`|x_tx_t^p|`$ as the probability of making an error. Note that the expectation is taken over the probabilistic prediction, whereas for the deterministic $`\mathrm{\Lambda }_\xi `$ algorithm the expectation is taken over the environmental distribution $`\mu `$. The multi-dimensional case \[KW99\] could then be interpreted as a (probabilistic) prediction of symbols over an alphabet $`𝒜=\{0,1\}^d`$, but error bounds for the absolute loss have yet to be proven. It would be interesting to generalize WM and bound (18) to arbitrary alphabet and to general loss functions with probabilistic interpretation.

### 5.5 Miscellaneous

Another direction is to investigate the learning aspect of universal prediction. Many prediction schemes explicitly learn and exploit a model of the environment. Learning and exploitation are melted together in the framework of universal Bayesian prediction. A separation of these two aspects in the spirit of hypothesis learning with MDL \[VL00\] could lead to new insights. The attempt at an information theoretic interpretation of Theorem 4 may be made more rigorous in this or another way. In the end, this may lead to a simpler proof of Theorem 4 and maybe even for the loss bounds. Finally, the system should be implemented and tested on specific induction problems for specific finite $`M`$ with computable $`\xi `$.

## 6 Summary

Solomonoff’s universal probability measure has been generalized to arbitrary probability classes and weights. A wise choice of $`M`$ widens the applicability by reducing the computational burden for $`\xi `$. A framework, where predictions result in losses of arbitrary, but known form, has been considered. Loss bounds for general loss functions have been proved, which show that the universal prediction scheme $`\mathrm{\Lambda }_\xi `$ can compete with the best possible informed scheme $`\mathrm{\Lambda }_\mu `$. The results show that universal prediction is ideally suited for games of chance with a sequence of bets, observations, and rewards. Extensions in various directions have been suggested.

## Acknowledgements

I want to thank Ray Solomonoff for many valuable discussions and for encouraging me to derive the general loss bounds presented here.


text:
# Learning multilayer perceptrons efficiently

## Abstract

A learning algorithm for multilayer perceptrons is presented which is based on finding the principal components of a correlation matrix computed from the example inputs and their target outputs. For large networks our procedure needs far fewer examples to achieve good generalization than traditional on-line algorithms.

PACS: 84.35.+i,89.80.+h, 64.60.Cn

Multilayer neural networks achieve the ability to approximate any reasonable function arbitrarily well  by interconnecting sufficiently many architecturally identical processing elements, the perceptrons. This replication of identical elements invariably leads to symmetries in the multilayer architecture which need to be broken during training to achieve good generalization . In the dynamics of gradient based training procedures the symmetries give rise to suboptimal fixed points and slow convergence. In particular this holds for stochastic gradient descent methods which to date have been the method of choice for training large networks because, at each time step, the update of the network parameters is based on the presentation of only a single training example and hence computationally cheap. Such stochastic gradient descent methods have been intensively studied in the framework of on-line learning, where each of the training examples is used just once by the network . This is attractive because there is no need to store the entire set of training examples, and the approach also works for nonstationary problems. However, in the on-line framework the number of required training examples is coupled to the slow temporal evolution caused by the symmetries and thus unreasonably large training sets are necessary. In fact the ratio of the number of examples needed for good generalization to the number of free parameters in the network diverges with the network size . While there have been investigations into optimizing the on-line dynamics , these have not lead to practical algorithms since the optimized procedures assume that the symmetry breaking provided by the initial conditions is macroscopic and known to the learning algorithm.

In this letter we present a learning algorithm which has many of the attractive features of the traditional on-line procedures but yields good generalization using a much smaller number of training examples. Further an exact analysis of the algorithm’s performance shows that the ratio of required examples to the number of free parameters in the networks does stay finite in the thermodynamic limit.

The multilayer architecture we analyze is a committee machine with $`K`$ hidden units defined by

$$\tau (\xi )=g\left(K^{1/2}\underset{i=1}{\overset{K}{}}h(B_i^T\xi )\right),$$
(1)
where $`\xi IR^N`$ is the input and the $`B_iIR^N`$ are the unknown parameter vectors. The goal of learning is to estimate these parameter vectors, which we shall refer to as teacher vectors, from a training set of $`P`$ examples $`(\xi ^\mu ,\tau (\xi ^\mu ))`$ of the input/output relationship. We shall initially focus on regression problems and later sketch the modifications to our procedure needed for classification tasks. For regression, the output function $`g`$ is usually assumed invertible, and this easily reduces to the linear case by applying the inverse function to the target outputs. So in this case we simply assume $`g(x)=x`$. For brevity we also assume that the $`B_i`$ are orthonormal.

In its simplest form our procedure can be seen as generalizing Hebbian learning. There the parameter vector $`\overline{J}^P`$ of a simple perceptron approximating the target function $`\tau (\xi )`$ is obtained as the average

$$\overline{J}^P=P^1\underset{\mu =1}{\overset{P}{}}\tau (\xi ^\mu )\xi ^\mu .$$
(2)
Our basic observation is that when $`\tau (\xi )`$ is given by a multilayer perceptron it is important to not only consider the empirical mean of the distribution of $`\tau (\xi ^\mu )\xi ^\mu `$ as in plain Hebbian learning but also its correlation matrix. In particular some of the eigenvectors of the correlation matrix correspond to the parameter vectors $`B_i`$ and thus the supervised learning task is mapped onto the well understood problem of principal component analysis. While this mapping by itself does not solve the original problem, it does provide a crucial reduction in its dimensionality, such that the remaining problem becomes almost trivial in the thermodynamic limit $`N\mathrm{}`$.

We shall consider the general correlation matrix

$$C^P=P^1\underset{\mu =1}{\overset{P}{}}F\left(\tau (\xi ^\mu )\right)\xi ^\mu \xi _{}^{\mu }{}_{}{}^{T}$$
(3)
where the simple choice $`F(x)=x^2`$ for the weight function $`F`$ yields the analogy to Hebbian learning. Assuming the components of the input vectors $`\xi ^\mu `$ to be independent Gaussian random variables with zero mean and unit variance, it is straightforward to analyze the spectrum of $`C^P`$ in the limit $`P\mathrm{}`$. The spectrum then only has three eigenvalues $`\lambda _0,\overline{\lambda }`$ and $`\lambda _\mathrm{\Delta }`$. The degeneracy of $`\lambda _0`$ is $`NK`$ and any vector orthogonal to all teacher vectors $`B_i`$ is an eigenvector. The eigenvalue $`\overline{\lambda }`$ has the single eigenvector $`\overline{B}=K^{1/2}_{i=1}^KB_i`$; this eigenvector is of little interest since for large $`P`$ one also has $`\overline{B}\overline{J}^P`$, and it is thus simpler to use Hebb’s rule (2) to estimate $`\overline{B}`$. The important eigenspace is the one of $`\lambda _\mathrm{\Delta }`$ since it is spanned by the $`K1`$ vectors $`B_1B_j`$ ($`j=2,\mathrm{},K`$). The three eigenvalues can be written as averages over independent, zero mean, unit variance, Gaussian random variables $`y_i`$. For instance one obtains $`\lambda _\mathrm{\Delta }=\frac{1}{2}F(\tau _y)(y_1y_2)^2_y`$ where $`\tau _y=K^{1/2}_{i=1}^Kh(y_i)`$. The activation function $`h`$ is sigmoidal (odd, monotonic and bounded), and when stating specific numerical values we will always assume $`h(y)=\text{erf}(y)`$. Then for the choice $`F(x)=x^2`$ one finds the ordering $`\overline{\lambda }>\lambda _0>\lambda _\mathrm{\Delta }`$ and $`\lambda _0\lambda _\mathrm{\Delta }=\frac{1}{K}\frac{8}{3\pi }(1\frac{1}{\sqrt{5}})`$.

For a finite number $`P`$ of training examples the degeneracy in the spectrum is broken by random fluctuations. But a computation of the orthonormal eigenvectors $`\mathrm{\Delta }_j^P`$ corresponding to the $`K1`$ smallest eigenvalues of $`C^P`$ nevertheless yields an estimate of the space spanned by the difference vectors $`B_1B_j`$. To measure the success of approximating this space, we introduce the overlap

$$\rho =(K1)^{1/2}\text{Tr}(\mathrm{\Delta }^TBB^T\mathrm{\Delta })^{1/2},$$
(4)
where B is the matrix $`(B_1,\mathrm{},B_K)`$ of the teacher vectors and $`\mathrm{\Delta }=(\mathrm{\Delta }_1^P,\mathrm{},\mathrm{\Delta }_{K1}^P)`$. This is a sensible measure because $`\rho `$ is invariant with respect to the choice of an orthonormal basis of the space spanned by the $`\mathrm{\Delta }_j^P`$ and since it attains its maximal value of $`1`$ iff all $`\mathrm{\Delta }_j^P`$ lie in the space spanned by the $`B_i`$. Simulations showing the overlap $`\rho `$ as function of the number of examples per free parameter, $`\alpha =P/(KN)`$, are depicted in Fig. 1. They indicate a second order phase transition from zero to positive $`\rho `$ at a critical value $`\alpha _\mathrm{c}`$. The evolution of $`\rho `$ can be calculated exactly in the thermodynamic limit $`N\mathrm{}`$. But since up to now we are only estimating the space spanned by the $`B_i`$, instead of the vectors themselves, we address this problem first and defer the calculation of $`\rho `$.

Using the above eigenspace procedure, the $`B_i`$ can be approximated by linear combinations of the $`\mathrm{\Delta }_j^P`$ and $`\overline{J}^P`$. Thus the original $`KN`$ dimensional problem is reduced to $`K^2`$ dimensions and this is much easier for large $`N`$. To exploit the dimensionality reduction, we write the possible linear combinations as $`\stackrel{~}{\mathrm{\Delta }}\mathrm{\Gamma }`$ where $`\stackrel{~}{\mathrm{\Delta }}=(\overline{J}^P,\mathrm{\Delta }_1^P,\mathrm{},\mathrm{\Delta }_{K1}^P)`$, $`\mathrm{\Gamma }`$ is a $`K`$ by $`K`$ matrix of tunable parameters, and we define a student network $`\sigma _\mathrm{\Gamma }`$ in the restricted space via $`\sigma _\mathrm{\Gamma }(\xi )=g\left(K^{1/2}_{i=1}^Kh\left((\stackrel{~}{\mathrm{\Delta }}\mathrm{\Gamma })_i^T\xi \right)\right)`$. We want to choose $`\mathrm{\Gamma }`$ to minimize the generalization error $`ϵ_\mathrm{g}=\frac{1}{2}(\tau (\xi )\sigma _\mathrm{\Gamma }(\xi ))^2_\xi `$ and to this end $`\widehat{P}`$ additional examples $`(\widehat{\xi }^\nu ,\tau (\widehat{\xi }^\nu ))`$, $`\nu =1,\mathrm{},\widehat{P}`$, are used. To simplify the theoretical analysis, the additional examples should be picked independently of the examples used to obtain $`\stackrel{~}{\mathrm{\Delta }}`$. We then apply the standard on-line gradient descent procedure $`\mathrm{\Gamma }^{\nu +1}=\mathrm{\Gamma }^\nu \eta _\mathrm{\Gamma }(\tau (\widehat{\xi }^\nu )\sigma _{\mathrm{\Gamma }^\nu }(\widehat{\xi }^\nu ))^2`$. However, by choosing a scaling of the learning rate $`\eta `$ such that $`\eta 1`$ for large $`N`$, the stochasticity drops out of the procedure  and in the restricted space it performs gradient descent in $`ϵ_\mathrm{g}`$. Further we can scale $`\widehat{P}`$ such that $`\eta \widehat{P}1`$ and then $`\mathrm{\Gamma }^{\widehat{P}}`$ will be a minimum of $`ϵ_\mathrm{g}`$ for large $`N`$. Finally both scaling conditions can be satisfied while observing $`\widehat{P}N`$, so that the required number of additional examples is negligible compared to the size of the first training set. Note that thanks to the reduced dimensionality the details of the scaling of $`\eta `$ and $`\widehat{P}`$ only affect finite size effects and not the performance in the thermodynamic limit. Simulations combining the two stages of our algorithm (Fig. 1) show a steady decay in $`ϵ_\mathrm{g}`$.

We now turn to the theoretical calculation of $`\rho `$ in the important first stage of the algorithm. The smallest eigenvalue of $`C^P`$ can be found by minimizing $`J^TC^PJ`$ under the constraint $`|J|=1`$. Hence we consider the partition function $`Z=𝑑J\mathrm{exp}(\beta PJ^TC^PJ)`$ where the integration is over the $`N`$-dimensional unit sphere. For large $`N`$ the typical properties of the minimization problem are found by calculating the training set average $`\mathrm{ln}Z`$ and taking the limit $`\beta \mathrm{}`$. Using replicas, one introduces the $`K`$ dimensional order parameter $`R`$ , the typical overlap $`B^TJ`$ of the teacher vectors with a vector $`J`$ drawn from the Gibbs distribution $`Z^1\mathrm{exp}(\beta PJ^TC^PJ)`$, and further the replica symmetric scalar order parameter $`q`$ which is the overlap $`J_{}^{1}{}_{}{}^{T}J^2`$ of two vectors drawn from this distribution. For $`P>N`$ the correlation matrix $`C^P`$ is nonsingular, $`q`$ approaches $`1`$ with increasing $`\beta `$, and the relevant scaling is $`\chi =\beta (1q)=𝒪(1)`$. Then for the quenched average of $`\mathrm{ln}Z`$ one finds in the limit $`\beta \mathrm{}`$

$$\frac{\mathrm{ln}Z}{\beta N}=\begin{array}{c}\text{extr}\\ R,\chi \end{array}R^TA(\alpha ,\chi )R+a(\alpha ,\chi ).$$
(5)
Here $`a(\alpha ,\chi )=\alpha KG_\chi (\tau _y)_y+\frac{1}{2\chi }`$, the $`K`$ by $`K`$ matrix $`A(\alpha ,\chi )`$ has entries

$$A_{jk}(\alpha ,\chi )=\alpha KG_\chi (\tau _y)(y_jy_k\delta _{jk})_y\frac{\delta _{jk}}{2\chi },$$
(6)
and $`G_\chi (\tau _y)=F(\tau _y)/(1+2\chi F(\tau _y))`$. Since Eq. (5) is quadratic in $`R`$ the extremal problem can only have a solution $`R0`$ if the matrix $`A`$ is singular. From the symmetries one easily obtains that $`A`$ has just two eigenvalues. The first can be written as $`A_{11}A_{12}`$, it is the relevant eigenvalue in our case  and its eigenspace is spanned by the vectors $`e_1e_j`$ ($`j=2,\mathrm{},K`$), where $`e_1,\mathrm{},e_K`$ denote the standard basis of $`IR^K`$. This degeneracy shows that the difference between the $`K1`$ smallest eigenvalues of the correlation matrix vanishes for large $`N`$. So the simple procedure of analyzing the properties of the single vector $`J`$ minimizing $`J^TC^PJ`$, in fact yields the properties of the $`K1`$ eigenvectors vectors of $`C_P`$ with smallest eigenvalues in the thermodynamic limit.

Due to the degeneracy, we can reparametrize (5) setting $`R=\rho (e_1e_j)/\sqrt{2}`$ and obtain an extremal problem with only two variables $`\rho `$ and $`\chi `$. Note that $`\rho `$ is indeed the parameter introduced in the analysis of the numerical simulations. Its evolution is now obtained by solving (5) and this confirms the continuous phase transition found in the simulations from $`\rho =0`$ to positive $`\rho `$ at a critical value $`\alpha _\mathrm{c}`$. For $`K=2`$ one finds $`\alpha _\mathrm{c}=4.49`$ and for $`K=3`$ the result is $`\alpha _\mathrm{c}=8.70`$. As shown in Fig. 1 beyond the phase transition there is excellent agreement between the theoretical prediction and the simulations.

To obtain generic results for large $`K`$ note that the contributions of $`y_1`$ and $`y_2`$ to the target output $`\tau _y`$ will be small in this limit. Decomposing $`\tau _y`$ as $`\tau _y=\tau _y^{}+\delta _y/\sqrt{K}`$ where $`\delta _y=h(y_1)+h(y_2)`$, and expanding $`G_\chi (\tau _y)`$ up to second order for large $`K`$ simplifies Eq. 5 to:

$`{\displaystyle \frac{\mathrm{ln}Z}{\beta N}}=\begin{array}{c}\text{extr}\\ \rho ,\chi \end{array}`$ $`{\displaystyle \frac{\alpha \rho ^2}{4}}G_\chi ^{\prime \prime }(\tau _y^{})_y\delta _y^2((y_1y_2)^22)_y`$ (10)
$`\alpha KG_\chi (\tau _y)_y+{\displaystyle \frac{1\rho ^2}{2\chi }}.`$
On the one hand, applying the central limit theorem, the multiple integrals $`G_\chi ^{\prime \prime }(\tau _y^{})_y`$ and $`G_\chi (\tau _y)_y`$ can now be replaced by a single average, on the other hand the structure of (10) shows that $`\chi `$ approaches zero with increasing $`K`$. These observations yield the very simple result that $`\rho ^2=1\alpha _\mathrm{c}/\alpha `$ for large $`K`$ and $`\alpha >\alpha _\mathrm{c}`$. The value of $`\alpha _\mathrm{c}`$ is obtained as

$$\alpha _\mathrm{c}=\frac{4KF^2(\mu z)_z}{F^{\prime \prime }(\mu z)_z^2\left(z^2h^2(z)_z\mu ^22zh(z)_z^2\right)^2}$$
(11)
where $`\mu ^2=h^2(z)_z`$, and the distribution of $`z`$ is Gaussian with zero mean and unit variance. It is now straightforward to derive the optimal choice of $`F`$ from Eq. (11) by noting that in the denominator $`F^{\prime \prime }(\mu z)_z=\mu ^2(z^21)F(\mu z)_z`$. Applying the Cauchy-Schwarz inequality to $`F^2(\mu z)_z/F^{\prime \prime }(\mu z)_z`$ then yields that the optimal choice is $`F(x)=x^2\mu ^2`$. For this choice the eigenvalue $`\lambda _0`$ of the correlation matrix $`C^P`$ vanishes for large $`P`$. So the optimal $`F`$ maximizes the signal to noise ratio between the eigenspace of difference vectors we want to estimate and the orthogonal space.

For the activation function $`h(x)=\text{erf}(x)`$ one finds that $`\alpha _\mathrm{c}=1.96K`$ when $`F`$ is optimal, whereas the simple choice $`F(x)=x^2`$ yields a critical $`\alpha `$ which is one and a half times higher. Simulation results for $`K=7`$ plotted in Fig. 2 show that the large $`K`$ theory provides a reasonable approximation already for networks with quite few hidden units.

To round off the analysis of the regression problem, we obtain a theoretical prediction for the generalization error achieved by combining the two stages of our procedure. A simple calculation shows that the overlap $`r`$ of the Hebbian vector $`\overline{J}^\mu `$ with $`\overline{B}`$, $`r=\overline{B}^T\overline{J}^\mu /|\overline{J}^\mu |`$, for large $`N`$ satisfies $`r=(1+\frac{\mathrm{arcsin}(2/3)}{2\alpha K})^{1/2}`$. Further, using the explicit expression for $`ϵ_\mathrm{g}`$ given in  and the values $`r`$ and $`\rho `$, we can calculate the minimum of $`ϵ_\mathrm{g}`$ in the restricted space and find a theoretical prediction for the generalization error obtained by the second stage. This prediction is compared to the simulations in Fig. 1 and 2.

We next consider classification problems, that is we assume that the output function of the network given by Eq. (1) is $`g(x)=\text{sgn}(x)`$. Then, since the output is binary, $`\lambda _0=\overline{\lambda }`$ holds for any choice of the weight function $`F`$, and our procedure cannot immediately be applied. However, it is possible to gain information about the target rule by not just considering its output $`\tau ^\mu `$ but by comparing the output to the value $`\overline{B}^T\xi ^\mu `$ which would have been obtained by the linearized network, $`g(x)=h(x)=x`$. This is feasible since, even for classification, $`\overline{B}`$ can be estimated by the Hebbian vector $`\overline{J}^P`$ defined by (2). We are thus lead to consider more general correlation matrices of the form $`C^P=P^1_{\mu =1}^PF(\tau (\xi ^\mu ),\overline{J}_{}^{\mu 1}{}_{}{}^{T}\xi ^\mu )\xi ^\mu \xi _{}^{\mu }{}_{}{}^{T}.`$ A reasonable way of choosing $`F`$ is to focus on inputs $`\xi ^\mu `$ where the target output has a different sign than the output of the linearized network. So we consider $`F(x,y)=\mathrm{\Theta }(xy)\mu `$, where $`\mathrm{\Theta }`$ is the Heavyside step function.

In the large $`P`$ limit, the matrix $`C^P`$ has the same eigenspaces as in the case of regression and the three eigenvalues will in general be different. For the activation function we chose $`h(x)=\text{sgn}(x)`$ to compare our results with the findings for on-line Gibbs learning , to our knowledge the only algorithm which has been simulated in any detail for classifications task with a connected committee machine. For this $`h`$ one finds $`\lambda _\mathrm{\Delta }>\lambda _0>\overline{\lambda }`$, and for large $`K`$ to leading order $`\lambda _\mathrm{\Delta }\lambda _0=\sqrt{2\pi 4}/(\pi ^2K)`$. Numerical simulations of the procedure are shown in Fig. 3 using $`\mu =\pi ^1\mathrm{arccos}\sqrt{2/\pi }`$. Motivated by our findings in the case of regression, this choice of $`\mu `$ yields $`\lambda _0=0`$ for large $`K`$. While the training sets needed to achieve good generalization are much larger than for regression, they are, for exactly the same architecture, approximately $`20`$-times smaller than for on-line Gibbs learning already for $`N=150`$.

Throughout this Letter we have assumed that the inputs are drawn from an isotropic distribution. Typically, in practical applications, the input data itself will have some structure, and in this case the inputs have to be whitened before applying our procedure. To this end one computes the correlations of the inputs by $`D^P=P^1_{\mu =1}^P\xi ^\mu \xi _{}^{\mu }{}_{}{}^{T}`$. Then instead of just considering $`C^P`$, one finds the extremal eigenvectors $`\stackrel{~}{\mathrm{\Delta }}_j`$ of the transformed matrix $`R_{}^{1}{}_{}{}^{T}C^PR^1`$, where the square matrix $`R`$ is the Cholesky factorization of $`D^P`$, $`R^TR=D^P`$. Then an estimate for the space spanned by the difference vectors is obtained by transforming back, setting $`\mathrm{\Delta }_j=R^1\stackrel{~}{\mathrm{\Delta }}_j`$. Numerical simulations for $`F(x)=x^2`$ show that whitening noticeably improves the performance of our algorithm even when the inputs are picked from an isotropic distribution. This is due to the fact that the spectrum of the empirical correlation matrix $`D^P`$ has a finite width even when it converges to a delta peak with increasing $`P`$.

In summary, we have presented the first learning algorithm for a realistic multilayer network which can be exactly analyzed in the thermodynamic limit and yields good generalization already for a finite ratio of training examples to free parameters. This contrasts with the behavior of traditional on-line procedures  for target functions such as the ones considered in this Letter. As long as there are on the order of $`N`$ on-line steps, the dynamics is exactly described by deterministic differential equations for a set of self-averaging order parameters  in the limit of large $`N`$. However, even when tuned by choosing a good learning rate, the dynamics is stuck in a badly generalizing plateau on this time scale unless some information about the target function is built into the initial conditions. For the realistic case of randomly chosen initial values, the symmetries are eventually broken due to small initial fluctuations and the accumulation of fluctuations during the dynamics, but only when the number of on-line steps is much larger than $`N`$. This symmetry breaking is not adequately described by the deterministic differential equations. When the training set is sampled without replacement, this divergence of the time scale means that large numbers of examples are needed. But even in the recently analyzed scenario of sampling with replacement , the above theoretical problems remain and are compounded since the analysis requires involved techniques such as dynamical replica theory and quite a few approximations already on the time scale $`N`$.

Hence, we believe that the findings in this Letter open new avenues for research and that in particular the performance of the algorithm for classification tasks deserves a more detailed analysis. But the perhaps most intriguing aspect of our procedure is that it does not really assume that the number of hidden units in the architecture is known a priori. This is highlighted by the inset of Fig. 2 showing that the number of hidden units can be determined by inspecting the eigenvalue spectrum of the correlation matrix. So our findings also provide a novel perspective on the problem of model selection.

The work of two of us (C.B. and R.U.) was supported by the Deutsche Forschungsgemeinschaft.