Paper deep dive
A mean teacher algorithm for unlearning of language models
Yegor Klochkov
Intelligence
Status: succeeded | Model: google/gemini-3.1-flash-lite-preview | Prompt: intel-v1 | Confidence: 96%
Last extracted: 3/11/2026, 1:10:30 AM
Summary
The paper introduces a mean teacher algorithm for language model unlearning, demonstrating that it approximates a trajectory of slow natural gradient descent (NGD). To address vanishing gradients associated with traditional unlearning losses, the author proposes a new loss function called 'negative log-unlikelihood' (NLUL). The combination of the mean teacher algorithm and NLUL is shown to improve performance on the MUSE benchmark, particularly in balancing memorization reduction with model utility preservation.
Entities (5)
Relation Signals (3)
Mean Teacher Algorithm → usedfor → Machine Unlearning
confidence 100% · A mean teacher algorithm for unlearning of language models
Negative Log-Unlikelihood → improvesperformanceon → MUSE Benchmark
confidence 95% · We show that the combination of mean teacher and NLUL improves some metrics on the MUSE benchmarks
Mean Teacher Algorithm → approximates → Natural Gradient Descent
confidence 90% · We show that the mean teacher can approximate a trajectory of a slow natural gradient descent (NGD)
Cypher Suggestions (0)
No Cypher suggestions yet.
Abstract
Abstract:One of the goals of language model unlearning is to reduce memorization of selected text instances while retaining the model's general abilities. Despite various proposed methods, reducing memorization of large datasets without noticeable degradation in model utility remains challenging. In this paper, we investigate the mean teacher algorithm (Tarvainen & Valpola, 2017), a simple proximal optimization method from continual learning literature that gradually modifies the teacher model. We show that the mean teacher can approximate a trajectory of a slow natural gradient descent (NGD), which inherently seeks low-curvature updates that are less likely to degrade the model utility. While slow NGD can suffer from vanishing gradients, we introduce a new unlearning loss called "negative log-unlikelihood" (NLUL) that avoids this problem. We show that the combination of mean teacher and NLUL improves some metrics on the MUSE benchmarks (Shi et al., 2024).
Tags
Links
Trouble viewing inline? Open PDF directly →
Full Text
184,337 characters extracted from source content.
Expand or collapse full text
A mean teacher algorithm for unlearning of language models Yegor Klochkov Abstract One of the goals of language model unlearning is to reduce memorization of selected text instances while retaining the model’s general abilities. Despite various proposed methods, reducing memorization of large datasets without noticeable degradation in model utility remains challenging. In this paper, we investigate the mean teacher algorithm (Tarvainen & Valpola, 2017), a simple proximal optimization method from continual learning literature that gradually modifies the teacher model. We show that the mean teacher can approximate a trajectory of a slow natural gradient descent (NGD), which inherently seeks low-curvature updates that are less likely to degrade the model utility. While slow NGD can suffer from vanishing gradients, we introduce a new unlearning loss called “negative log-unlikelihood” (NLUL) that avoids this problem. We show that the combination of mean teacher and NLUL improves some metrics on the MUSE benchmarks (Shi et al., 2024). Machine Learning, ICML 1 Introduction In the recent years, Large Language Models (LLMs) have reached unprecedented capabilities that are achieved through training on massive datasets. Often comprising hundreds of billions of tokens, these datasets are collected from diverse internet sources. This data-intensive approach raises significant ethical concerns as some protected user and copyright data may be compromised. The General Data Protection Regulation (GDPR) establishes fundamental rights for individuals, including the “right to be forgotten,” which grants the deletion of personal data upon request. Concurrently, the inclusion of copyrighted content in training datasets has initiated legal quarrels (Doe 1 v. GitHub, 2022; Tremblay v. OpenAI, 2023). These regulatory and legal actions have sparked interest in methods that would reduce the influence of selective training data on predictions of large pretrained models, a process known as machine unlearning (Bourtoule et al., 2021; Ginart et al., 2019). In unlearning of language models, we are aimed at modifying the model weights that would make them behave as if it was never trained on particular selected data, which we refer to as forget set (Yao et al., 2023; Eldan & Russinovich, 2023; Cooper et al., 2024; Shi et al., 2024). For example, we want to suppress certain outputs, like being able to reproduce some memorized training instances word-by-word with completion requests, or answering questions about their content. Following (Shi et al., 2024), we refer to the completion ability as verbatim memorization, and the ability to answer questions as knowledge memorization. We therefore want to reduce such memorization, ideally, without decreasing the model’s utility. There have been a lot of interest in the recent literature towards language model unlearning, with various methods proposed. A generally accepted approach consists of finetuning a model by optimizing an objective specifically designed to discourage memorization. Much of research is dedicated to designing specific loss functions that encourage unlearning, with notable examples such as log-likelihood (referred to as gradient ascent) (Thudi et al., 2022; Yao et al., 2023), Negative Preference Optimization (Zhang et al., 2024a), and others (Chundawat et al., 2023; Fan et al., 2024; Wang et al., 2024; Li et al., 2024). With the resulting objective including a regularization term, it is typically optimized with AdamW (Loshchilov et al., 2017), which is a standard way of optimizing language models. A few recent works dedicated to unlearning address the question of how we optimize the corresponding objective, Wang et al. (2025) employ some gradient projections, and Bu et al. (2024) consider adaptive learning rates designed specifically for unlearning. We revisit a simple method called mean teacher (Tarvainen & Valpola, 2017), which so far appeared only in the context of continual learning literature. We show that in certain scenarios, it can approximate the trajectory of a slow natural gradient descent. This means that the algorithm performs updates along the low curvature directions, which can help retain the model utility and avoid neural collapse. However, we show that mean teacher can suffer from vanishing gradients with some popular choices of unlearning losses, such as log-likelihood and NPO. Instead, we introduce a very simple loss called Negative Log UnLikelihood (NLUL) that avoids this problem and combines well with mean teacher. For the experiments, we focus on the MUSE benchmark (Shi et al., 2024), which contains rather large forget set (3.3M tokens). Although there are methods that are capable of removing verbatim memorization, the reported results so far show it is challenging to reduce knowledge memorization without dramatic reduction in the utility. We show that one of the versions of mean teacher can achieve this. However, we show that in this case the reduction of knowledge memorization is accompanied by reduction of the Massive Multitask Language Understanding metric (MMLU) (Hendrycks et al., 2020). This means that the evaluations in the MUSE-benchmark may not be enough for adequate assessment of unlearning algorithms that researchers state reduce knowledge memorization. On the other hand, there are variants that achieve competititve verbatim memorization reduction and utility preservation, while improving on another metrics that are part of the MUSE benchmark associated with risks of privacy leakage. In our work, we focus exclusively on methods that edit the model weights. We mention that some methods are based on expanding the model decision space through guardrails (Thaker et al., 2024b; Liu et al., 2024) and training low-rank adapters on top of the pretrained model weights (Gao et al., 2024; Ji et al., 2024). Unlearning is sometimes seen as a safety mechanism, with focus on reducing potentially harmful knowledge of language models (Li et al., 2024). There, the focus is on removing knowledge of whole topics or concepts. On the other hand, we focus on unlearning independent text instances. 2 Unlearning set-up and description of the method Most unlearning benchmarks assume that we are given an access to two datasets: • forget set DfsubscriptD_fDitalic_f is a set of training instances that we know are the source of knowledge we want to remove; • retain set DrsubscriptD_rDitalic_r is another set of training instances that do not contain undesirable knowledge, which is intended to assist in retaining the model’s utility. We acknowledge knowing the source of undesirable knowledge would already be a strong assumption. However, most of the existing benchmarks follow such set up (Maini et al., 2024; Shi et al., 2024), see also a review in (Thaker et al., 2024a). Often, unlearning is performed by optimizing an objective of the form, Lforget(Df;θ)+Lretain(Dr;θ),subscriptsubscriptsubscriptsubscriptL_forget(D_f;θ)+L_retain(D_r;θ),Litalic_f o r g e t ( Ditalic_f ; θ ) + Litalic_r e t a i n ( Ditalic_r ; θ ) , (1) where LforgetsubscriptL_forgetLitalic_f o r g e t averages some unlearning loss on the given examples. That is, Lforget(Df;θ)=1|Df|∑(x,y)∈Dfℓ(h(x;θ),y),subscriptsubscript1subscriptsubscriptsubscriptℓℎL_forget(D_f;θ)= 1|D_f| _(x,y)∈ D_f (h(x;θ% ),y),Litalic_f o r g e t ( Ditalic_f ; θ ) = divide start_ARG 1 end_ARG start_ARG | Ditalic_f | end_ARG ∑( x , y ) ∈ D start_POSTSUBSCRIPT f end_POSTSUBSCRIPT ℓ ( h ( x ; θ ) , y ) , (2) where h(x;θ)ℎh(x;θ)h ( x ; θ ) are the logit outputs of the model on the given instance x, and y is the target label. One popular choice of unlearning loss, is the log likelihood (L) ℓLL(h,y)=logsf(h)y,subscriptℓℎsfsubscriptℎ _L(h,y)= (h)_y,ℓitalic_L L ( h , y ) = log sf ( h )y , where sf(h)y=ehy/(∑jehj)sfsubscriptℎsuperscriptsubscriptℎsubscriptsuperscriptsubscriptℎsf(h)_y=e^h_y/ ( _je^h_j )sf ( h )y = eitalic_hitalic_y / ( ∑j eitalic_hitalic_j ) is the softmax function. Such loss is a straightforward negation of the standard learning loss—negative log likelihood (NLL)—and it is intended to discourage correct prediction on the inputs from the forget set. Other unlearning losses include Negative Preference Optimization (NPO) (Zhang et al., 2024a); The Incompetent Teacher loss outputs KL divergence between the outputs h(x;θ)ℎh(x;θ)h ( x ; θ ) and that of a model of smaller size that does not exhibit memorization (Chundawat et al., 2023). In Section 4, we additionally introduce a new loss called negative log-unlikelihood. While NPO is designed to address the problem with explosive gradients when unlearning with L, our new loss also addresses vanishing gradients at the start of the unlearning process. In Section 4, we elaborate on this issue in detail, and we also provide some empirical arguments in Section 5. The penalizing term LretainsubscriptL_retainLitalic_r e t a i n is designed to maintain the utility of the model, which is usually defined as either NLL loss on the instances from DrsubscriptD_rDitalic_r, or the Kullback-Leibler (KL) divergence with the outputs of the initial model. To be precise, in the latter case the retain loss looks as follows, Lretain(Dr;θ)=1|Dr|∑(x,y)∈DrKL(sf(h(x;θ))∥sf(h(x;θ0)))subscriptsubscript1subscriptsubscriptsubscriptconditionalsfℎsfℎsubscript0L_retain(D_r;θ)= 1|D_r| _(x,y)∈ D_rKL(sf(% h(x;θ))\;\|\;sf(h(x; _0)))Litalic_r e t a i n ( Ditalic_r ; θ ) = divide start_ARG 1 end_ARG start_ARG | Ditalic_r | end_ARG ∑( x , y ) ∈ D start_POSTSUBSCRIPT r end_POSTSUBSCRIPT K L ( sf ( h ( x ; θ ) ) ∥ sf ( h ( x ; θ0 ) ) ) (3) The experiments often suggest that KL regularization leads to better results (Shi et al., 2024; Yao et al., 2023). Everywhere in the paper we assume that the datasets contain next tokens as prediction labels, i.e., DfsubscriptD_fDitalic_f is formed of pairs (x,y)(x,y)( x , y ) — context and next token — where x=(s0,…,st−1)subscript0…subscript1x=(s_0,…,s_t-1)x = ( s0 , … , sitalic_t - 1 ), y=stsubscripty=s_ty = sitalic_t, and s is a pretraining sequence with t<|s|t<|s|t < | s |. Similarly, DrsubscriptD_rDitalic_r consists of context and next token pairs for sequences from the retain set. 3 Mean teacher algorithm and it’s approximation of natural gradient descent When the retain loss as in (3) is used, the optimization problem (1) can be seen as a proximal optimization problem minθL(θ)+(θ,θ′),subscriptsuperscript′ _θL(θ)+D(θ,θ ),minitalic_θ L ( θ ) + D ( θ , θ′ ) , where we have L(θ)=Lforget(θ)subscriptL(θ)=L_forget(θ)L ( θ ) = Litalic_f o r g e t ( θ ), the reference model θ′=θ0superscript′subscript0θ = _0θ′ = θ0 is equal to the starting point, and the divergence term is the KL loss on the retain set (θ,θ′)=1|Dr|∑(x,y)∈DrKL(sf(h(x;θ))∥sf(h(x;θ′))).superscript′1subscriptsubscriptsubscriptconditionalsfℎsfℎsuperscript′D(θ,θ )= 1|D_r| _(x,y)∈ D_rKL(% sf(h(x;θ))\;\|\;sf(h(x;θ ))).D ( θ , θ′ ) = divide start_ARG 1 end_ARG start_ARG | Ditalic_r | end_ARG ∑( x , y ) ∈ D start_POSTSUBSCRIPT r end_POSTSUBSCRIPT K L ( sf ( h ( x ; θ ) ) ∥ sf ( h ( x ; θ′ ) ) ) . (4) The divergence term measures the proximity of the model θ to the reference θ′θ θ′. For example, it is non-negative, and in fact, in the case of KL divergence it is also locally quadratic, in the sense that there is a matrix H(θ)H(θ)H ( θ ), such that (θ,θ′)=12(θ−θ′)⊤H(θ′)(θ−θ′)+o(‖θ−θ′‖2)superscript′12superscriptsuperscript′topsuperscript′superscriptnormsuperscript′2D(θ,θ )= 12(θ-θ ) % H(θ )(θ-θ )+o(\|θ-θ \|^2)D ( θ , θ′ ) = divide start_ARG 1 end_ARG start_ARG 2 end_ARG ( θ - θ′ )⊤ H ( θ′ ) ( θ - θ′ ) + o ( ∥ θ - θ′ ∥2 ) (5) In the case of KL divergence, this matrix corresponds to the Gauss-Newton Hessian (GNH), H(θ)=1|Dr|∑x∈Dr[Jθh(x;θ)]Sh(x;θ)[Jθh(x;θ)]⊤,1subscriptsubscriptsubscriptdelimited-[]subscriptℎsubscriptℎsuperscriptdelimited-[]subscriptℎtopH(θ)= 1|D_r| _x∈ D_r[J_θh(x;θ)]S_h(x;% θ)[J_θh(x;θ)] ,H ( θ ) = divide start_ARG 1 end_ARG start_ARG | Ditalic_r | end_ARG ∑x ∈ D start_POSTSUBSCRIPT r end_POSTSUBSCRIPT [ Jitalic_θ h ( x ; θ ) ] Sitalic_h ( x ; θ ) [ Jitalic_θ h ( x ; θ ) ]⊤ , (6) where Jθh(θ)=(∂hj(θ)∂θi)ijsubscriptℎsubscriptsubscriptℎsubscriptJ_θh(θ)= ( ∂ h_j(θ)∂ _i% )_ijJitalic_θ h ( θ ) = ( divide start_ARG ∂ hitalic_j ( θ ) end_ARG start_ARG ∂ θitalic_i end_ARG )i j is the Jacobian, and Sh=Diag(sf(h))−sf(h)sf(h)⊤subscriptℎDiagsfℎsfℎsfsuperscriptℎtopS_h=Diag(sf(h))-sf(h)sf(h) Sitalic_h = Diag ( sf ( h ) ) - sf ( h ) sf ( h )⊤. The Gauss-Newton Hessian is semi-positive definite by design, and it is often regarded as an approximation to the Hessian of the NLL loss (Martens, 2020). We also note that the formula above is only a particular case of the GNH corresponding to the cross entropy loss, and it can be defined for other losses as well. Furthermore, in this particular case, it is equivalent to the Fisher Information Matrix (Schraudolph, 2002; Kunstner et al., 2019). Let us now describe our version of mean teacher algorithm. For the sake of generality, let us assume that we have some loss function L(θ)L(θ)L ( θ ), that is not necessarily an unlearning loss, and some divergence (θ,θ′)superscript′D(θ,θ )D ( θ , θ′ ) that is not necessarily the KL divergence. The mean teacher works in gradient steps and maintains the model we are optimizing θtsubscript _tθitalic_t along with the reference model θt′subscript′ _t θitalic_t′ that is also updated with every iteration as exponentially weighted average of the past models, θt′=(κη)θt+(κη)(1−κη)θt−1+(κη)(1−κη)2θt−2+…,superscriptsubscript′subscript1subscript1superscript12subscript2… _t =(κη) _t+(κη)(1-κη) _t% -1+(κη)(1-κη)^2 _t-2+…,θitalic_t′ = ( κ η ) θitalic_t + ( κ η ) ( 1 - κ η ) θitalic_t - 1 + ( κ η ) ( 1 - κ η )2 θitalic_t - 2 + … , The two models are initialized with the given starting model θ0′=θ0superscriptsubscript0′subscript0 _0 = _0θ0′ = θ0, followed by the updates 1. θt+1←θt−η∇αL(θt)+(θt,θt′)←subscript1subscript∇subscriptsubscriptsuperscriptsubscript′ _t+1← _t-η∇\α L( _t)+D(% _t, _t )\θitalic_t + 1 ← θitalic_t - η ∇ α L ( θitalic_t ) + D ( θitalic_t , θitalic_t′ ) ; 2. θt+1′←(1−ηκ)θt′+ηκθt+1←superscriptsubscript1′1superscriptsubscript′subscript1 _t+1 ←(1-ηκ) _t +ηκ% _t+1θitalic_t + 1′ ← ( 1 - η κ ) θitalic_t′ + η κ θitalic_t + 1, where η is the learning rate, and κ is the contraction hyperparameter, and we additionally introduce a weight α<11α<1α < 1 for the loss term. The first step is a gradient descent step on the regularized loss, while the second updates the reference model by slightly sliding it towards the optimized model. Choosing a small weight α can ensure that the optimized model θtsubscript _tθitalic_t will stay close to the reference model θt′subscript′ _t θitalic_t′, while the sliding reference update ensures that we can still progress away from the initial point θ0subscript0 _0θ0. Suppose that the given divergence (θ,θ′)superscript′D(θ,θ )D ( θ , θ′ ) satisfies (5) with some H(θ)H(θ)H ( θ ). With some abuse of terms, we refer to this matrix as Hessian. In addition, we consider adding a quadratic regularization to improve the condition number of the Hessian λ(θ,θ′):=(θ,θ′)+λ2‖θ−θ′‖2,Hλ(θ):=H(θ)+λI.formulae-sequenceassignsubscriptsuperscript′2superscriptnormsuperscript′2assignsubscriptD_λ(θ,θ ):=D(θ,θ^% )+ λ2\|θ-θ \|^2, H_λ(% θ):=H(θ)+λ I.Ditalic_λ ( θ , θ′ ) := D ( θ , θ′ ) + divide start_ARG λ end_ARG start_ARG 2 end_ARG ∥ θ - θ′ ∥2 , Hitalic_λ ( θ ) := H ( θ ) + λ I . (7) Furthermore, we additionally apply momentum accumulation to the gradients (Polyak, 1964; Nemirovskij & Yudin, 1983). The pseudo code of the resulting algorithm is summarized in Algorithm 1. Below we show that for a sufficiently small coefficient α, the trajectory of Algorithm 1 follows a gradient descent, where the matrix H(θ)H(θ)H ( θ ) serves as a conditioner. In the case, where H(θ)H(θ)H ( θ ) is GNH, such trajectory corresponds to natural gradient descent (Amari, 1998). Algorithm 1 Mean teacher algorithm Input: Learning rate η, damping parameter λ, momentum μ, weight α, number of steps T, contraction parameter κ, starting parameter θ0subscript0 _0θ0 Initialize reference θ0′:=θ0assignsuperscriptsubscript0′subscript0 _0 := _0θ0′ := θ0. for t=11t=1t = 1 to T do θt←θt−1−η∇θαL(θt−1)+λ(θt−1,θt−1′)+μ(θt−1−θt−2)←subscriptsubscript1subscript∇subscript1subscriptsubscript1superscriptsubscript1′subscript1subscript2 _t← _t-1-η _θ \α L( _t-% 1)+D_λ( _t-1, _t-1 ) \+μ(% _t-1- _t-2)θitalic_t ← θitalic_t - 1 - η ∇θ α L ( θitalic_t - 1 ) + Ditalic_λ ( θitalic_t - 1 , θitalic_t - 1′ ) + μ ( θitalic_t - 1 - θitalic_t - 2 ) θt′←(1−ηκ)θt−1′+ηκθt←superscriptsubscript′1superscriptsubscript1′subscript _t ←(1-ηκ) _t-1 +ηκ% _tθitalic_t′ ← ( 1 - η κ ) θitalic_t - 1′ + η κ θitalic_t end for Return: θTsubscript _Tθitalic_T Theorem 3.1. Set γ:=καη/(1−κη)assign1γ:=καη/(1-κη)γ := κ α η / ( 1 - κ η ) and λ¯:=λ+(1−μ)κ/(1−ηκ)assign¯11 λ:=λ+(1-μ)κ/(1-ηκ)over¯ start_ARG λ end_ARG := λ + ( 1 - μ ) κ / ( 1 - η κ ). Consider the following updates θ¯t+1=θ¯t−γHλ¯(θ¯t)−1∇L(θ¯t−1),θ¯0=θ0formulae-sequencesubscript¯1subscript¯subscript¯superscriptsubscript¯1∇subscript¯1subscript¯0subscript0 θ_t+1= θ_t-γ H_ λ(% θ_t)^-1∇ L( θ_t-1)\,, % θ_0= _0over¯ start_ARG θ end_ARGt + 1 = over¯ start_ARG θ end_ARGt - γ Hover¯ start_ARG λ end_ARG ( over¯ start_ARG θ end_ARGt )- 1 ∇ L ( over¯ start_ARG θ end_ARGt - 1 ) , over¯ start_ARG θ end_ARG0 = θ0 (8) Suppose that 1) κ is a positive constant, 2) η, α are sufficiently small, 3) number of steps T is such that TγT γ is bounded by a constant, 4) H(θ)H(θ)H ( θ ) is symmetric positive-definite and satisfies Eq. (5) 5) (θ,θ′)superscript′D(θ,θ )D ( θ , θ′ ), L(θ)L(θ)L ( θ ), and H(θ)H(θ)H ( θ ) satisfy some regularity conditions, which we postpone to the appendix (Condition D.2). Then, the trajectory of Algorithm 1 approximately matches that of (8), with the following upper-bound maxt≤T‖θt−θ¯t‖≤O(αlog(1/α)).subscriptnormsubscriptsubscript¯1 _t≤ T\| _t- θ_t\|≤ O(α (1/α))\,.maxitalic_t ≤ T ∥ θitalic_t - over¯ start_ARG θ end_ARGt ∥ ≤ O ( α log ( 1 / α ) ) . Note, that the contraction coefficient κ effectively increases the damping parameter. This effect can be reduced with the use of momentum. For example, for a typical choice μ=0.90.9μ=0.9μ = 0.9, we get λ¯≈λ+0.1κ¯0.1 λ≈λ+0.1 ¯ start_ARG λ end_ARG ≈ λ + 0.1 κ. We additionally note that the condition on TγT γ bounds the total sum of step sizes in (8). The proof is deferred to Section D in the appendix. Alternative divergences The KL divergence is the most common way to regularize in proximal optimization methods. We additionally consider the Quadratic KL(QKL) loss defined as (θ,θ′)superscript′ (θ,θ )D ( θ , θ′ ) =1|Dr|∑(x,y)∈DrQKL(h(x;θ),h(x;θ′))absent1subscriptsubscriptsubscriptℎsuperscript′ = 1|D_r| _(x,y)∈ D_rQKL(h(x;θ),h(x;θ% ))= divide start_ARG 1 end_ARG start_ARG | Ditalic_r | end_ARG ∑( x , y ) ∈ D start_POSTSUBSCRIPT r end_POSTSUBSCRIPT Q K L ( h ( x ; θ ) , h ( x ; θ′ ) ) (9) QKL(h,h′)ℎsuperscriptℎ′ QKL(h,h )Q K L ( h , h′ ) :=(h−h′)⊤(Diag(sf(h))−sf(h)sf(h)⊤)(h−h′),assignabsentsuperscriptℎsuperscriptℎ′topDiagsfℎsfℎsfsuperscriptℎtopℎsuperscriptℎ′ :=(h-h ) (Diag(sf(h))-sf% (h)sf(h) )(h-h ),:= ( h - h′ )⊤ ( Diag ( sf ( h ) ) - sf ( h ) sf ( h )⊤ ) ( h - h′ ) , which attains the same approximation as the KL divergence in Eq. (5), with GNH in place of H(θ)H(θ)H ( θ ). We show in the experiments section that using QKL can sometimes give different results and improve some of the metrics in the benchmarks. For the sake of example, we also mention that H(θ)H(θ)H ( θ ) does not always have to be the GNH. If we have a loss function L(θ)L(θ)L ( θ ) that is convex, we can consider the associated Bregman divergence, (θ,θ′)=L(θ)−L(θ′)−(θ−θ′)⊤∇L(θ′),superscript′superscriptsuperscript′top∇superscript′D(θ,θ )=L(θ)-L(θ )-(θ-% θ ) ∇ L(θ ),D ( θ , θ′ ) = L ( θ ) - L ( θ′ ) - ( θ - θ′ )⊤ ∇ L ( θ′ ) , which due to Taylor expansion allows quadratic approximation with the original Hessian H(θ)=∇2L(θ).superscript∇2H(θ)=∇^2L(θ).H ( θ ) = ∇2 L ( θ ) . Such divergence was is in (Amid et al., 2022) for a private mirror descent algorithm. Using pretraining data for calculating divergence Tarvainen & Valpola (2017) introduce mean teacher as a semi-supervised finetuning method with the focus on classification models. The idea is to use a limited labeled set for the loss, and a large set of unlabeled examples for the divergence, which does not depend on the labels. In language model finetuning, all text examples are self-supervised since we always predict the next token. However, the examples we select still undergo some curation. For example, in MUSE-Books the forget contains Harry Potter books, while the retain split contains articles from Harry Potter Fandom Wiki pages. If we are faced with a new unlearning task, we will have to obtain curate a new retain split that matches the topic of a forget set in a similar way. Furthermore, in the unlearning benchmarks, the retain split is often not much larger and sometimes smaller than the forget set. In order to make the set up closer to the semi-supervised setting, we propose to replace the retain split with a subset of pretraining data. For language models, we propose to use the OpenWebText-2, which is readily available on Huggingface111https://huggingface.co/datasets/Skylion007/openwebtext. This dataset has moderate size by modern standards and often included in larger pretraining corpora (Gao et al., 2020). The dataset comprises of <<<1B tokens and our algorithms only use a fraction of that by sampling an independent batch with each gradient update. Apart from the advantage of having more tokens than a benchmark’s retain split, it can also be reused when we are faced with a new unlearning task which reduces the cost of collecting and curating the data. We additionally note that using pretraining subset DptsubscriptD_ptDitalic_p t instead of DrsubscriptD_rDitalic_r in Eq. (6) may be better aligned with preserving general abilities of the model. Ghorbani et al. (2019) numerically demonstrated that the subspaces spanned by top eigenvectors of GNH are aligned with that of Hessian of the cross entropy loss H(θ)=∇2LCE(θ)superscript∇2subscriptH(θ)=∇^2L_CE(θ)H ( θ ) = ∇2 Litalic_C E ( θ ). The natural gradients reduce the influence of these directions due to inversion of the Hessian, which therefore ensures that our updates have low curvature w.r.t. the pretraining loss LCE(θ;Dpt)=1|Dpt|∑(x,y)∈DptℓCE(h(x,θ),y).subscriptsubscript1subscriptsubscriptsubscriptsubscriptℓℎL_CE(θ;D_pt)= 1|D_pt| _(x,y)∈ D_pt _CE(h(x,% θ),y).Litalic_C E ( θ ; Ditalic_p t ) = divide start_ARG 1 end_ARG start_ARG | Ditalic_p t | end_ARG ∑( x , y ) ∈ D start_POSTSUBSCRIPT p t end_POSTSUBSCRIPT ℓitalic_C E ( h ( x , θ ) , y ) . This can help better preserve the general abilities of the language model that are not directly associated with the unlearning task at hand. We therefore replace DrsubscriptD_rDitalic_r with the pretraining subset DptsubscriptD_ptDitalic_p t in the experiments for two reasons: 1) it is a bigger set than DrsubscriptD_rDitalic_r and it does not depend on the unlearning task at hand 2) it can help us better preserve general abilities of the model. Implementation details Finally, notice that Theorem 3.1 assumes full-batch gradient updates. In the experiments, we make the following changes according to standard practice: 1. Batching: The gradients are calculated on batches of data, i.e., for each update we sample a batch of sequences from DfsubscriptD_fDitalic_f and a batch DptsubscriptD_ptDitalic_p t. The loss and the KL divergence are calculated over that batch. 2. Gradient clipping: we apply norm clipping to the gradients to avoid explosive gradients. Note, since that clipping is equivalent to reducing the learning rate value, we correspondingly reduce the contracting parameter κ←κ/max(‖g‖/c,1)←norm1κ←κ/ (\|g\|/c,1)κ ← κ / max ( ∥ g ∥ / c , 1 ), where c is the clipping coefficient. 3. Momentum: the gradients are fed into the momentum SGD optimizer. For reference, we provide the pseudo-code for this version of the algorithm in the appendix, Section A. 4 Unlearning losses Here we describe three popular unlearning losses, in addition to one which we introduce in this paper. The given objectives are to be minimized. 1. Log likelihood (L) is a straightforward negation of the standard learning loss — negative log likelihood (NLL), ℓLL(h;y)=logsf(h)y,subscriptℓℎsfsubscriptℎ _L(h;y)= (h)_y,ℓitalic_L L ( h ; y ) = log sf ( h )y , which is also sometimes referred to as gradient ascent, since it negates the gradients that would be used in learning with NLL. Using this loss for unlearning has two caveats: firstly, we typically start from a point where the we are pretty much converged on DforgetsubscriptD_forgetDitalic_f o r g e t due to memorization, and it is hard to escape this point, unless we are using an accelerated optimization method such as AdamW. However, this does not solve the second problem: we will observe explosive gradients the further away we get from the minima, which can eventually lead to a collapse of the model. 2. Negative preference optimization (NPO) is a modification of direct preference optimization loss (DPO), that only includes the negative examples, which are the ones sampled from the forget set, ℓNPOβ(s;θ)=−2βlogσ(−βlogπθ(s)πbase(s))subscriptℓsubscript2subscriptsubscript _NPO_β(s;θ)=- 2β σ (-β % _θ(s) _base(s) )ℓitalic_N P O start_POSTSUBSCRIPT β end_POSTSUBSCRIPT ( s ; θ ) = - divide start_ARG 2 end_ARG start_ARG β end_ARG log σ ( - β log divide start_ARG πitalic_θ ( s ) end_ARG start_ARG πitalic_b a s e ( s ) end_ARG ) (Zhang et al., 2024a) introduce this loss to address the problem with explosive gradients. In particular, they observe that the gradients of this loss correspond to reweighed gradients of the log-likelihood, 1|s|∇θℓNPOβ(s;θ)=π(s)βπθ(s)β+πbase(s)β1|s|∑(x,y)∈s∇θℓLL(h(x,θ);y).1subscript∇subscriptℓsubscriptsuperscriptsubscriptsuperscriptsubscriptsuperscript1subscriptsubscript∇subscriptℓℎ 1|s| _θ _NPO_β(s;θ)= π(s)^β% _θ(s)^β+ _base(s)^β 1|s| _(x,y)∈ s% _θ _L(h(x,θ);y)\,.divide start_ARG 1 end_ARG start_ARG | s | end_ARG ∇θ ℓitalic_N P O start_POSTSUBSCRIPT β end_POSTSUBSCRIPT ( s ; θ ) = divide start_ARG π ( s )β end_ARG start_ARG πitalic_θ ( s )β + πitalic_b a s e ( s )β end_ARG divide start_ARG 1 end_ARG start_ARG | s | end_ARG ∑( x , y ) ∈ s ∇θ ℓitalic_L L ( h ( x , θ ) ; y ) . Here, the weight will approach zero, as the probabilities πθ(s)subscript _θ(s)πitalic_θ ( s ) on the forget set reduce. The term πbase(s)subscript _base(s)πitalic_b a s e ( s ) is the probability of the forget token according to the starting model, which we denote as “base” to avoid confusion with the reference model that we update according to the mean teacher rule. This means that we have to additionally store and forward propagate through one more model, which moderately increases the cost of unlearning. To avoid confusion, we remark that this loss does not conform with the formula Eq. 2 as the the weight is calculated per sequence. 3. Incompetent teacher (IT) is a simple and intriguing approach that was introduced in (Chundawat et al., 2023), but appears have been overlooked in the LLM unlearning literature. They propose to minimize the KL divergence between outputs of the optimized model and that of a fixed smaller model that potentially does not exhibit memorization — the incompetent teacher. For example, if we unlearn a LLama-2 7B model, we can “push” its outputs towards that of the TinyLlama 1.1B on the forget data, rather than maximizing the NLL. The loss in that case looks as follows, ℓIT(h;y)=KL(sf(h)∥sf(hit))subscriptℓℎconditionalsfℎsfsubscriptℎ _IT(h;y)=KL (sf(h)\;\|\;sf(h_it) )ℓitalic_I T ( h ; y ) = K L ( sf ( h ) ∥ sf ( hitalic_i t ) ) where hitsubscriptℎh_ithitalic_i t are the logit outputs of the incompetent teacher. KL divergence, being convex w.r.t. the logits hℎh, is better suited to combine with popular accelerated optimization techniques, such as AdamW. Notice that the mismatch loss from (Yao et al., 2023) is equivalent to incompetent teacher choosing uniformly random tokens. 4. Negative log-unlikelihood (NLUL) is another candidate that we introduce in this paper, and it looks as follows ℓNLUL(h;y)=−log(1−sf(h)y)subscriptℓℎ1sfsubscriptℎ _NLUL(h;y)=- (1-sf(h)_y)ℓitalic_N L U L ( h ; y ) = - log ( 1 - sf ( h )y ) The idea is similar to the L approach, but instead of minimizing the likelihood, we maximize the log-unlikelihood log(1−p(y))1 (1-p(y))log ( 1 - p ( y ) ). Similar to NPO, such reformulation allows a better reweighing of the gradients per token ∇hℓNLUL(h;y)=sf(h)y1−sf(h)y∇hℓLL(h;y).subscript∇ℎsubscriptℓℎsfsubscriptℎ1sfsubscriptℎsubscript∇ℎsubscriptℓℎ _h _NLUL(h;y)= sf(h)_y1-sf(h)_y% _h _L(h;y)\,.∇h ℓitalic_N L U L ( h ; y ) = divide start_ARG sf ( h )y end_ARG start_ARG 1 - sf ( h )y end_ARG ∇h ℓitalic_L L ( h ; y ) . It has a similar effect when p(y)≈00p(y)≈ 0p ( y ) ≈ 0, helping to avoid explosive gradients. However, it additionally allows to escape the stagnating starting point by making more aggressive updates in the beginning, when we often have p(y)≈11p(y)≈ 1p ( y ) ≈ 1 in the forget set. It does not introduce additional hyperparameters (such as β in NPO) and does not require to maintain an additional base model. In the experiments, we find this loss combines best with MT. 5 Experiments 5.1 Benchmarks We use two benchmarks from (Shi et al., 2024), MUSE-News and MUSE-Books. In both cases, the authors provide us with a forget set, a retain set, a set of queries and ground truth labels for testing memorization and utility, and a target model to be unlearned. MUSE-News is based on BBC news articles collected after the release date of Llama-2. These articles are split into forget and retrain sets randomly. Then, they finetune a Llama-2 7B on both sets, so that it memorizes both sets of articles. The goal is to reduce memorization of the forget data while retaining knowledge of the retain articles. For this, they provide a set of completion and question queries extracted from the news articles, which provides the measurement of both retain and forget articles. Overall, the following metrics are considered: • verbmem_f: this metric measures verbatim memorization of the forget data – the ability to reproduce part of the article when prompted with the completion request. For this, they provide a 100 pairs prompt-completion, and measure the ROUGE-L score between the response to a prompt and ground truth completion. • knowmem_f: this metric evaluates knowledge memorization of the forget data– the ability to answer question about some facts contained in the article. For this, Shi et al. (2024) provide 100 pairs of questions and ground truth answers extracted from the forget split articles. Then, we measure how similar model’s response to the ground truth answer. They propose to use ROUGE-L to measure this similarity. • knowmem_r: similarly to knowledge memorization of the forget data, this metric measures that of the retain data, with QA pairs extracted from the retain articles. This value is proposed as a measure of utility of the model. • PrivLeak: they additionally propose a metric of how easy it is to detect that the forget set was used for unlearning with some state-of-the-art membership inference methods for language models (Shi et al., 2023). For exact definition, see Section 3.1 in (Shi et al., 2024), and we use the code they provide for evaluations. This metric can have both positive and negative values, and in the ideal situation it is close to zero. In addition to that, we measure the MMLU accuracy (Hendrycks et al., 2020), which is a popular metric of language model’s general knowledge. While original set of questions is composed of approximately 150K questions, to speed up the evaluations, we use a small validation set consisting of ≈1500absent1500≈ 1500≈ 1500 questions. We denote the corresponding accuracy as mmlu_val. In the MUSE-Books benchmark, the forget set contains four Harry Potter books, while the retain set consists of random articles from the Harry Potter FanWiki page. They similarly construct a set of queries, allowing one to measure verbmem_f, knowmem_f, knowmem_r, and we additionally include mmlu_val. The target model that was finetuned on both forget and retain splits, and it is a version of Llama-2 architecture as well. In both cases, we discard the provided retain split and instead rely on the pretraining dataset OpenWebText. 5.2 Methods For baselines, we optimize unlearning loss penalized by KL divergence on DptsubscriptD_ptDitalic_p t with AdamW optimizer. The methods are referred to as AdamW + (L/NPO/IT/NLUL) + KL. The proposed methods utilize mean teacher, with one of the four unlearning losses, and either KL or QKL regularization on DptsubscriptD_ptDitalic_p t. The mean teacher algorithm is implemented in accordance with Algorithm 2 in the appendix. The corresponding methods are labeled as MT + (L/NPO/IT/NLUL) + (KL/QKL). For the IT loss, we use TinyLlama-1.1B222https://huggingface.co/TinyLlama/TinyLlama-1.1B-Chat-v1.0, which shares the tokenizer with Llama-2. We list all hyperparameters for each experiment in Section B in the appendix. 5.3 Mean teacher can suffer from vanishing gradients Figure 1: NLL loss on the forget set for MT using different unlearning losses (green). For IT we only show the KL divergence between the bad teacher and the target model. For NPO/L we additionally perform 2 epochs with AdamW to “escape” the starting point. Figure 2: Gradient norms during MT training in Figure 2. We first perform experiments with each of the four losses to find out which combines best with MT. We use MUSE-News data and report the training loss in Figure 2. In the cases of L, NPO, NLUL, we report the NLL loss on the forget set. In the case of IT we report KL divergence with the logit outputs of TinyLlama-1.1B. Observe that in the latter case the algorithm does not manage to reduce memorization even after 3000 steps, despite the fact that the KL loss keeps reducing. Furthermore, when used with NPO and L losses, the method struggles to escape the initial point, see the green line on the bottom plots in Figure 2. To show that this is only a problem in the starting point, we perform a 2 epoch finetuning with AdamW, and once we escaped that starting point with low loss, MT manages to continue increasing NLL. On the contrary, NLUL manages to increase NLL without any additional tricks. We show the gradient norm per step in Figure 2, which also shows lower order gradients magnitudes for L and NPO than NLUL. We conclude that despite a plausible natural gradient approximation, MT only approximates slow natural gradient descent and can suffer from vanishing gradients. Out of these four unlearning losses, NLUL combines best with MT out of these four choices, and we will use the MT + NLUL combination in the remaining experiments. Table 1: Results for mean teacher and baselines on MUSE-News/Books. We run each method 3 times and report average and std of each metric, and we also report number of epochs. The (⋆ ⋆) sign indicates that the method is run until knowmem_f≤31.1absent31.1≤ 31.1≤ 31.1. PrivLeak is only measured on one run. verbmem_f ↓ ↓ knowmem_f ↓ ↓ knowmem_r ↑ ↑ PrivLeak mmlu_val ↑ ↑ MUSE-News Target 56.3 64.2 55.2 -99.8 40.0 Retrain 21.8 33.1 55.0 0.0 NA MT + NLUL + KL 14.8 (2.3) 58.2 (0.5) 52.7 (0.9) -38.8 39.0 (0.5) MT + NLUL + QKL 8.4 (0.7) 40.0 (6.2) 50.8 (0.5) -6.9 31.0 (4.3) AdamW + NPO + KL 0.0 (0.0) 51.1 (1.3) 50.1 (0.9) 89.0 39.1 (0.8) AdamW + IT + KL 21.3 (0.2) 64.4 (0.1) 52.7 (0.2) -99.5 40.5 (0.1) AdamW + L + KL 0.0 (0.0) 51.2 (1.2) 50.7 (0.5) 91.9 38.0 (0.0) AdamW + NLUL + KL 2.9 (0.7) 58.4 (0.3) 49.9 (0.9) 94.2 39.6 (0.2) MT + NLUL + QKL⋆\,\! ⋆ 8.5 (0.9) 24.8 (3.1) 51.3 (0.2) 25.0 19.5 (0.6) MUSE-Books Target 99.8 59.4 66.9 -57.5 26.3 Retrain 14.3 28.9 74.5 0.0 NA NLUL + MT + KL 1.5 (0.3) 22.0 (3.8) 59.1 (0.7) -35.9 26.6 (0.2) NLUL + MT + QKL 8.8 (1.0) 17.5 (1.5) 43.0 (3.2) -57.7 25.6 (0.1) AdamW + NPO + KL 0.0 (0.0) 17.7 (6.8) 37.17 (5.3) -24.1 25.8 (0.4) Table 2: Examples of answers from target model of the MUSE-News benchmark and it’s unlearned modifications. We show one example probing verbatim memorization (completion) and one example probing knowledge memorization (QA). Examples are cherry-picked for illustration purpose. Prompt / Question […] The BMA wants a 35% pay rise to make up for 15 years of below-inflation wage rises, It says falling pay has caused a recruitment and retention crisis. Dr Vivek Trivedi, co-director of the junior doctors’ committee at the BMA, said on Saturday that Mr Barclay is yet to put a serious offer on the table. ”All we’re asking for is a credible offer that shows us he’s serious, that we can start a path of negotiations to Which three nuclear power plants were taken offline in Germany by midnight on Saturday? GT try to address the real-terms pay cut,” he said. The strike is due to take place from 07:00 BST on Tuesday until 07:00 BST on Saturday. The NHS national medical director, Professor Sir Stephen Powis, is warning that the strike will cause ”unparalleled levels of disruption” as it is longer than previous strikes and comes after the bank holiday when many staff are ”taking much-needed holiday”. Up to quarter of a million operations and appointments could be postponed because Isar 2, Emsland and Neckarwestheim 2 MT + NLUL + KL hopefully everybody getting everybody what they want,” he said Unterscheidung between junior doctors and everybody else is ”not helping anybody”, he added. The strike is due to take place from 07:00 BST on Tuesday until 07:00 on Saturday. nobody will be nobody will be nobody will be nobody will be nobody will be nobody will be nobody will be nobody will be nobody will be nobody will be […] 3 - Isar 2, Emsland and Neckarwestheim 2 MT + NLUL + QKL∗\,\!^*∗ try to address everybody’s concerns,” he said Unterscheidung zwischen den verschiedenen Arten von Arzneimitteln, die in der Apotheke erhältlich sind, und deren Anwendung. ”If we don’t get that, then there’s no point asking junior doctors to go back to work because they’l just be walking out again in a few weeks time.” The BMA says junior doctors represent everybody from those fresh out of medical school to nobody with more than everybody in between. The union says it has the support of more than 70% of the medical profession kwietynka, Is sierpnik, Neckarwestheim nobody knows how many reactors there are at each plant AdamW + NPO + KL ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” ” 3 - Isar 2, Emsland and Neckarwestheim 2 Figure 3: Comparison of mean teacher (blue), baselines with pretraining data (red), original baselines from (Shi et al., 2024) that use retain split for regularization (green), and recent new methods as reported in (Bu et al., 2024; Wang et al., 2024; Fan et al., 2024; Wang et al., 2025) (orange). 5.4 Comparison to baselines Let us first discuss the results on the MUSE-News benchmark. We compare the performance of MT with baselines in Table 1. For baselines, we include four options AdamW + (NPO/IT/L/NLUL) + KL. Following (Shi et al., 2024), we run the algorithms for a few epochs until we reach verbatim memorization ≤21.8absent21.8≤ 21.8≤ 21.8, which corresponds to the retrained model’s value. Most of the methods show moderate reduction in the utility (knowmem_r). In Figure 3 we additionally show results reported in (Shi et al., 2024) for the baselines that use the retain split (red), along with some of the more recent methods333As of beginning of April 2025, out of 45 citations of (Shi et al., 2024) we only found 4 papers that: a) introduce a new method b) test it on MUSE-News.: NGDiff (Bu et al., 2024), GRU (Wang et al., 2025), FLAT (Wang et al., 2024), and SimNPO (Fan et al., 2024) (orange). The baselines from (Shi et al., 2024) show greater utility loss, see leftmost graph in Figure 3, which suggests that using the pretraining set helps retain the utility better. Next, we look at how the methods balance utility and reducing knowledge memorization on the forget set (knowmem_f), see the middle plot in Figure 3. One of the points that (Shi et al., 2024) make is that there is no free lunch when it comes to reducing knowledge memorization (see their Figure 5), and it always comes at almost proportional reduction in utility (knowmem_r). On the first glance, our results suggest that the method MT + NLUL + QKL makes such attempt. Furthermore, we run the unlearning algorithm a little longer, until we reach the measurement knowmem_f≤33.1absent33.1≤ 33.1≤ 33.1, which corresponds to the model that was only finetuned on the retain split and never seen any of the forget data. This result is reported as MT + NLUL + QKL⋆ in Table 1 and Figure 3. We can see that knowmem_f is further reduced down to 24.8, while knowmem_r stays close the level of most competitive methods. However, as reported in Table 1, for both of this methods we observe significant reduction in MMLU evaluation (mmlu_val). Therefore, the no free lunch observation of (Shi et al., 2024) still stands. We additionally show in Section 5.5 below, that further finetuning on completely unrealted data can revert both knowmem_f and mmlu_val back. This means that the knowledge is still somehow encoded in the weights after the unlearning is performed. These observations suggests that researchers who report reduction in knowmem_f should report similar evaluations. Finally, we report that MT-based algorithm produce models with lowered privacy leakage, see column PrivLeak in Table 1 (recall that for this metric, closer to zero is better). This is also true when we compare with more recent works, as shown in the rightmost plot in Figure 3. We note that NGDiff and GRU (Bu et al., 2024; Wang et al., 2025) do not report PrivLeak. In Section C in the appendix, we additionally report the the MUSE-News sustainability experiment. We run it on MT + NLUL + QKL and AdamW + NPO + KL, and both demonstrate good ability to perform sequential unlearning requests. Finally, we highlight that on MUSE-Books, both Mean Teacher implementations reduce verbatim and knowledge memorization without significant harm to both knowmem_r and MMLU accuracy. However, we note the concerning factor in this case, that the MMLU performance is sufficiently degraded already in the target model. Table 3: Recovering knowledge memorization on MUSE-News after unlearning is performed. Results for unlearning with NLUL + MT + QKL∗\,\!^*∗, followed by finetuning with SFT and Magpie-Pro-300K-Filtered. unlearning SFT verbmem_f ↓ ↓ 12.3 14.3 knowmem_f ↓ ↓ 40.2 58.1 knowmem_r ↑ ↑ 51.6 53.8 mmlu_val ↑ ↑ 33.2 40.8 5.5 Discussion of knowledge memorization unlearning The output of MT + NLUL + QKL∗\,\!^*∗ on question answering prompts (see example in Table 2) and the reduced mmlu_val evaluation may give us a hint that the model general ability of answering questions has reduced. In order to restore the answering abilities, we finetune the unlearned model on the SFT dataset Magpie-Pro-300K-Filtered444https://huggingface.co/datasets/Magpie-Align/Magpie-Pro-300K-Filtered. Surprisingly, it recovers both MMLU accuracy and knowledge memorization (Table 3). Since the intersection between the alignment data and news articles is unlikely, we can argue that the knowledge still remained in the model weights after the unlearning procedure. Such occurrences are not uncommon, with some reporting restoring memorization with quantization (Zhang et al., 2024b), finetunning on a fraction of the forget set (Hu et al., 2024), and by intervening on internal representations (MacLaurin, 2024). 5.6 Conclusion on the experiments We sum up our observations as follows: 1. Mean teacher is on par with the state-of-the-art method in terms of the trade-off between verbatim memorization and utility on MUSE benchmarks. 2. Baselines’ performances improve when the pretraining data is utilized for regularization. 3. Although mean teacher is capable to reduce the knowledge memorization on the forget data without dramatic drop in knowledge or retain articles, additional evaluations show reduction of MMLU. Furthermore, we show that finetuning on unrelated alignment data can revert both of these metrics. We hope that in the future, a similar assessment will be made by researchers stating knowledge removal on MUSE-News. 4. Mean teacher algorithms can reduce privacy leakage compared to the baselines and results reported in the other recent papers (Shi et al., 2024; Wang et al., 2024; Fan et al., 2024). Overall, we have demonstrated that the proposed mean teacher algorithm is a competitive unlearning algorithm that attains low verbatim memorization, has good utility preservation, and reduces the risks of privacy leakage. 5.7 Reproducibility To facilitate reproducibility of our results, we have made our implementation available at https://github.com/yklochkov-bytedance/mt-unlearn . Acknowledgments We thank Zafar Takhirov for encouragement and helpful discussions. References Agarwal et al. (2017) Agarwal, N., Bullins, B., and Hazan, E. Second-order stochastic optimization for machine learning in linear time. Journal of Machine Learning Research, 18(116):1–40, 2017. Amari (1998) Amari, S.-I. Natural gradient works efficiently in learning. Neural computation, 10(2):251–276, 1998. Amid et al. (2022) Amid, E., Ganesh, A., Mathews, R., Ramaswamy, S., Song, S., Steinke, T., Suriyakumar, V. M., Thakkar, O., and Thakurta, A. Public data-assisted mirror descent for private model training. In International Conference on Machine Learning, p. 517–535. PMLR, 2022. Bourtoule et al. (2021) Bourtoule, L., Chandrasekaran, V., Choquette-Choo, C. A., Jia, H., Travers, A., Zhang, B., Lie, D., and Papernot, N. Machine unlearning. In 2021 IEEE Symposium on Security and Privacy (SP), p. 141–159. IEEE, 2021. Bu et al. (2024) Bu, Z., Jin, X., Vinzamuri, B., Ramakrishna, A., Chang, K.-W., Cevher, V., and Hong, M. Unlearning as multi-task optimization: A normalized gradient difference approach with an adaptive learning rate. arXiv preprint arXiv:2410.22086, 2024. Chundawat et al. (2023) Chundawat, V. S., Tarun, A. K., Mandal, M., and Kankanhalli, M. Can bad teaching induce forgetting? unlearning in deep networks using an incompetent teacher. In Proceedings of the AAAI Conference on Artificial Intelligence, volume 37, p. 7210–7217, 2023. Cooper et al. (2024) Cooper, A. F., Choquette-Choo, C. A., Bogen, M., Jagielski, M., Filippova, K., Liu, K. Z., Chouldechova, A., Hayes, J., Huang, Y., Mireshghallah, N., et al. Machine unlearning doesn’t do what you think: Lessons for generative ai policy, research, and practice. arXiv preprint arXiv:2412.06966, 2024. Doe 1 v. GitHub (2022) Doe 1 v. GitHub, I. 4:22-cv-06823. N.D. Cal., 2022. Eldan & Russinovich (2023) Eldan, R. and Russinovich, M. Who’s harry potter? approximate unlearning in llms. arXiv preprint arXiv:2310.02238, 2023. Fan et al. (2024) Fan, C., Liu, J., Lin, L., Jia, J., Zhang, R., Mei, S., and Liu, S. Simplicity prevails: Rethinking negative preference optimization for llm unlearning. arXiv preprint arXiv:2410.07163, 2024. Fisher et al. (2023) Fisher, J., Liu, L., Pillutla, K., Choi, Y., and Harchaoui, Z. Influence diagnostics under self-concordance. In International Conference on Artificial Intelligence and Statistics, p. 10028–10076. PMLR, 2023. Gao et al. (2024) Gao, C., Wang, L., Weng, C., Wang, X., and Zhu, Q. Practical unlearning for large language models. arXiv preprint arXiv:2407.10223, 2024. Gao et al. (2020) Gao, L., Biderman, S., Black, S., Golding, L., Hoppe, T., Foster, C., Phang, J., He, H., Thite, A., Nabeshima, N., et al. The pile: An 800gb dataset of diverse text for language modeling. arXiv preprint arXiv:2101.00027, 2020. Ghorbani et al. (2019) Ghorbani, B., Krishnan, S., and Xiao, Y. An investigation into neural net optimization via hessian eigenvalue density. In International Conference on Machine Learning, p. 2232–2241. PMLR, 2019. Ginart et al. (2019) Ginart, A., Guan, M., Valiant, G., and Zou, J. Y. Making ai forget you: Data deletion in machine learning. Advances in neural information processing systems, 32, 2019. Hendrycks et al. (2020) Hendrycks, D., Burns, C., Basart, S., Zou, A., Mazeika, M., Song, D., and Steinhardt, J. Measuring massive multitask language understanding. arXiv preprint arXiv:2009.03300, 2020. Hu et al. (2024) Hu, S., Fu, Y., Wu, S., and Smith, V. Jogging the memory of unlearned models through targeted relearning attacks. In ICML 2024 Workshop on Foundation Models in the Wild, 2024. Ji et al. (2024) Ji, J., Liu, Y., Zhang, Y., Liu, G., Kompella, R. R., Liu, S., and Chang, S. Reversing the forget-retain objectives: An efficient llm unlearning framework from logit difference. arXiv preprint arXiv:2406.08607, 2024. Klochkov & Liu (2024) Klochkov, Y. and Liu, Y. Revisiting inverse hessian vector products for calculating influence functions. arXiv preprint arXiv:2409.17357, 2024. Koh & Liang (2017) Koh, P. W. and Liang, P. Understanding black-box predictions via influence functions. In International conference on machine learning, p. 1885–1894. PMLR, 2017. Kunstner et al. (2019) Kunstner, F., Hennig, P., and Balles, L. Limitations of the empirical fisher approximation for natural gradient descent. Advances in neural information processing systems, 32, 2019. Li et al. (2024) Li, N., Pan, A., Gopal, A., Yue, S., Berrios, D., Gatti, A., Li, J. D., Dombrowski, A.-K., Goel, S., Phan, L., et al. The wmdp benchmark: Measuring and reducing malicious use with unlearning. arXiv preprint arXiv:2403.03218, 2024. Liu et al. (2024) Liu, C. Y., Wang, Y., Flanigan, J., and Liu, Y. Large language model unlearning via embedding-corrupted prompts. arXiv preprint arXiv:2406.07933, 2024. Loshchilov et al. (2017) Loshchilov, I., Hutter, F., et al. Fixing weight decay regularization in adam. arXiv preprint arXiv:1711.05101, 5:5, 2017. MacLaurin (2024) MacLaurin, D. Unlearning via RMU is mostly shallow. https://w.lesswrong.com/posts/6QYpXEscd8GuE7BgW/unlearning-via-rmu-is-mostly-shallow, 2024. Accessed: 2024-01-31. Maini et al. (2024) Maini, P., Feng, Z., Schwarzschild, A., Lipton, Z. C., and Kolter, J. Z. Tofu: A task of fictitious unlearning for llms. arXiv preprint arXiv:2401.06121, 2024. Martens (2020) Martens, J. New insights and perspectives on the natural gradient method. Journal of Machine Learning Research, 21(146):1–76, 2020. Nemirovskij & Yudin (1983) Nemirovskij, A. S. and Yudin, D. B. Problem complexity and method efficiency in optimization. 1983. Polyak (1964) Polyak, B. T. Some methods of speeding up the convergence of iteration methods. Ussr computational mathematics and mathematical physics, 4(5):1–17, 1964. Schraudolph (2002) Schraudolph, N. N. Fast curvature matrix-vector products for second-order gradient descent. Neural computation, 14(7):1723–1738, 2002. Shi et al. (2023) Shi, W., Ajith, A., Xia, M., Huang, Y., Liu, D., Blevins, T., Chen, D., and Zettlemoyer, L. Detecting pretraining data from large language models. arXiv preprint arXiv:2310.16789, 2023. Shi et al. (2024) Shi, W., Lee, J., Huang, Y., Malladi, S., Zhao, J., Holtzman, A., Liu, D., Zettlemoyer, L., Smith, N. A., and Zhang, C. MUSE: Machine unlearning six-way evaluation for language models. arXiv preprint arXiv:2407.06460, 2024. Tarvainen & Valpola (2017) Tarvainen, A. and Valpola, H. Mean teachers are better role models: Weight-averaged consistency targets improve semi-supervised deep learning results. Advances in neural information processing systems, 30, 2017. Thaker et al. (2024a) Thaker, P., Hu, S., Kale, N., Maurya, Y., Wu, Z. S., and Smith, V. Position: Llm unlearning benchmarks are weak measures of progress. arXiv preprint arXiv:2410.02879, 2024a. Thaker et al. (2024b) Thaker, P., Maurya, Y., and Smith, V. Guardrail baselines for unlearning in llms. arXiv preprint arXiv:2403.03329, 2024b. Thudi et al. (2022) Thudi, A., Deza, G., Chandrasekaran, V., and Papernot, N. Unrolling sgd: Understanding factors influencing machine unlearning. In 2022 IEEE 7th European Symposium on Security and Privacy (EuroS&P), p. 303–319. IEEE, 2022. Tremblay v. OpenAI (2023) Tremblay v. OpenAI, I. 23-cv-03416-amo. N.D. Cal., 2023. Wang et al. (2024) Wang, Y., Wei, J., Liu, C. Y., Pang, J., Liu, Q., Shah, A. P., Bao, Y., Liu, Y., and Wei, W. Llm unlearning via loss adjustment with only forget data. arXiv preprint arXiv:2410.11143, 2024. Wang et al. (2025) Wang, Y., Wang, Q., Liu, F., Huang, W., Du, Y., Du, X., and Han, B. Gru: Mitigating the trade-off between unlearning and retention for large language models. arXiv preprint arXiv:2503.09117, 2025. Yao et al. (2023) Yao, Y., Xu, X., and Liu, Y. Large language model unlearning. arXiv preprint arXiv:2310.10683, 2023. Zhang et al. (2019) Zhang, G., Li, L., Nado, Z., Martens, J., Sachdeva, S., Dahl, G., Shallue, C., and Grosse, R. B. Which algorithmic choices matter at which batch sizes? insights from a noisy quadratic model. Advances in neural information processing systems, 32, 2019. Zhang et al. (2024a) Zhang, R., Lin, L., Bai, Y., and Mei, S. Negative preference optimization: From catastrophic collapse to effective unlearning. arXiv preprint arXiv:2404.05868, 2024a. Zhang et al. (2024b) Zhang, Z., Wang, F., Li, X., Wu, Z., Tang, X., Liu, H., He, Q., Yin, W., and Wang, S. Catastrophic failure of LLM unlearning via quantization. arXiv preprint arXiv:2410.16454, 2024b. Appendix A Batched algorithm with momentum See pseudo-code in Algorithm 2. Note that (θ,θ′)superscript′D(θ,θ )D ( θ , θ′ ) in the pseudo-code can either denote KL, Eq. (4), or QKL, Eq. (9). We also consider introducing additional damping parameter λ according to Eq. 7. Algorithm 2 Batched mean teacher with norm clipping and momentum Input: Learning rate η, weight α, number of steps T, contraction parameter κ, starting parameter θ0subscript0 _0θ0, clip value c, batch size Initialize reference θ0′:=θ0assignsuperscriptsubscript0′subscript0 _0 := _0θ0′ := θ0 Initialize momentum v=00v=0v = 0 for t=11t=1t = 1 to T do Sample batch BfsubscriptB_fBitalic_f from DfsubscriptD_fDitalic_f randomly Sample batch BptsubscriptB_ptBitalic_p t from DptsubscriptD_ptDitalic_p t randomly gt←∇θαL(θt−1;Bf)+(θt−1,θt−1′;Bpt)←subscriptsubscript∇subscript1subscriptsubscript1superscriptsubscript1′subscriptg_t← _θ\α L( _t-1;B_f)+D(% _t-1, _t-1 ;B_pt)\gitalic_t ← ∇θ α L ( θitalic_t - 1 ; Bitalic_f ) + D ( θitalic_t - 1 , θitalic_t - 1′ ; Bitalic_p t ) Scaling factor after clipping l=1/max(‖gt‖,c)1normsubscriptl=1/ (\|g_t\|,c)l = 1 / max ( ∥ gitalic_t ∥ , c ) v←μv+lg←v←μ v+lgv ← μ v + l g θt←θt−1−ηv←subscriptsubscript1 _t← _t-1-η vθitalic_t ← θitalic_t - 1 - η v θt′←(1−lηκ)θt−1′+lηκθt←superscriptsubscript′1superscriptsubscript1′subscript _t ←(1-lηκ) _t-1 +lηκ% _tθitalic_t′ ← ( 1 - l η κ ) θitalic_t - 1′ + l η κ θitalic_t end for Return: θTsubscript _Tθitalic_T Appendix B Hyperparameters We use the following set of hyperparameters: • MT + (NPO/IT/NLUL/L) + (KL/QKL): η=0.00050.0005η=0.0005η = 0.0005, α=0.050.05α=0.05α = 0.05 for MUSE-News and α=0.10.1α=0.1α = 0.1 for MUSE-Books, λ=0.50.5λ=0.5λ = 0.5, κ=10.010.0κ=10.0κ = 10.0, batch size 40404040; • AdamW + (NPO/IT/NLUL/L) + KL: learning rate 0.000010.000010.000010.00001, β=(0.9,0.95)0.90.95β=(0.9,0.95)β = ( 0.9 , 0.95 ), warmup schedule with 100100100100 steps at 10%percent1010\%10 % + 100100100100 steps increasing linearly, α=0.050.05α=0.05α = 0.05 for MUSE-News and α=0.10.1α=0.1α = 0.1 for MUSE-Books, batch size 40404040; In each case we take one epoch as 100 steps, which with batch size 40404040 approximately equals one pass over MUSE-News. In each case the model is trained until __ verbmem\_ftypewriter_verbmem _ typewriter_f is lower or equal to the value of the retrained model. With the exception of MT+NLUL+KL⋆, where we train until __ knowmem\_ftypewriter_knowmem _ typewriter_f becomes lower than retrained model’s value. In Table 4 we report the number of epochs required for this. We report average and std over 3 runs. Table 4: Number of epochs for each method to reach required metrics verbmem_f/knowmem_f # epochs MUSE-News MT + NLUL + KL 4 MT + NLUL + QKL 4 AdamW + NPO + KL 2 AdamW + IT + KL 7.3 (0.3) AdamW + L + KL 3 AdamW + NLUL + KL 3 MT + NLUL + QKL⋆\,\! ⋆ 5.7 (1.0) MUSE-Books NLUL + MT + KL 15 NLUL + MT + QKL 13 AdamW + NPO + KL 3 Appendix C Sustainability experiment on MUSE-News Shi et al. (2024) propose and additional experiment to test how unlearning algorithms respond to sequential unlearning request which may happen the such algorithms are utilized in the industrial setting. For this, the forget split of the MUSE-News benchmark is split into forget1, forget2, forget3, forget4, each comprised into approximately 0.8M tokens. The baselines that they consider show poor utility preservation when faced with such sequential unlearning requests, see Figure 6b in (Shi et al., 2024). Here we conduct this experiment for the algorithms MT+NLUL+QKL and AdamW + NPO +KL, both using pretraining split for regularization, unlike the NPO implementation in (Shi et al., 2024). The original paper does not clarify the stopping rule for these sequential requests. In our attempt to conduct a fair experiment, we propose to match the number of steps for each algorithm as reported in Table 4. Since the datasets are smaller, the resulting number of epochs is in fact larger for each little split. We do so for both of the algorithm. We report the results in Figure 4. This shows that both the baseline and MT unlearning algorithm preserve utility well. This further highlights the importance of using a bigger pretraining set for regularization. Figure 4: Sustainability of unlearning: how utility preserves with sequential unlearning requests. We perform experiment for mean teacher and NPO, both using pretraining data instead of MUSE-News retain split. Appendix D Proof of Theorem 3.1 The proof of Theorem 3.1 is based on comparison to the inverse Hessian-vector product (IHVP) estimation algorithm (Agarwal et al., 2017; Koh & Liang, 2017). For example, in the case μ=00μ=0μ = 0, we show that the sequence approximately satisfies the iterative updates, ut+1≈ut−ηα∇L(θt)+Hλ¯(θt)ut+μ(ut−ut−1)subscript1subscript∇subscriptsubscript¯subscriptsubscriptsubscriptsubscript1u_t+1≈ u_t-η\α∇ L( _t)+H_ λ(% _t)u_t\+μ(u_t-u_t-1)uitalic_t + 1 ≈ uitalic_t - η α ∇ L ( θitalic_t ) + Hover¯ start_ARG λ end_ARG ( θitalic_t ) uitalic_t + μ ( uitalic_t - uitalic_t - 1 ) (10) For fixed θtsubscript _tθitalic_t, these updates would correspond to ones used for estimating IHVP, with the Hessian being Hκ¯(θt)subscript¯subscriptH_ κ( _t)Hover¯ start_ARG κ end_ARG ( θitalic_t ), and the vector being α∇L(θt)∇subscriptα∇ L( _t)α ∇ L ( θitalic_t ). By carefully taking into account the changes in θtsubscript _tθitalic_t, we can show that starting from some early iteration, θt−θt′=−αHλ¯(θt)−1∇L(θt)+O~(α2),subscriptsuperscriptsubscript′subscript¯superscriptsubscript1∇subscript~superscript2 _t- _t =-α H_ λ( _t)^-1% ∇ L( _t)+ O(α^2),θitalic_t - θitalic_t′ = - α Hover¯ start_ARG λ end_ARG ( θitalic_t )- 1 ∇ L ( θitalic_t ) + over~ start_ARG O end_ARG ( α2 ) , where O~(⋅)~⋅ O(·)over~ start_ARG O end_ARG ( ⋅ ) hides some logarithmic factors. This approximation then allows us to compare directly to the natural gradient descent trajectory. D.1 Momentum lemma for quadratic loss We first describe a convergence bound for inverse Hessian vector product (IHVP) algorithm that accounts for the presence of momentum and potential error, e.g. as in the approximate equation (10). IHVP algorithms can be described as gradient descent iterations applied to the quadratic objective (Agarwal et al., 2017; Koh & Liang, 2017; Fisher et al., 2023; Klochkov & Liu, 2024) minu12u⊤(H+λ)u+g⊤u,subscript12superscripttopsuperscripttop _u 12u (H+λ)u+g u,minitalic_u divide start_ARG 1 end_ARG start_ARG 2 end_ARG u⊤ ( H + λ ) u + g⊤ u , which attains minimum at u⋆=−(H+λ)−1gsuperscript⋆superscript1u =-(H+λ)^-1gu⋆ = - ( H + λ )- 1 g. Below we show a convergence bound for erroneous gradient descent with momentum. We note that similar results appear, e.g., in (Zhang et al., 2019). Lemma D.1. Suppose, we have a sequence of updates ut+1subscript1 u_t+1uitalic_t + 1 =ut−(ηg+η(H+λ)ut+ηϵt)+μ(ut−ut−1),t≥1,u−1:=u0formulae-sequenceabsentsubscriptsubscriptsubscriptitalic-ϵsubscriptsubscript1formulae-sequence1assignsubscript1subscript0 =u_t-(η g+η(H+λ)u_t+η _t)+μ(u_t% -u_t-1), t≥ 1, u_-1:=u_0= uitalic_t - ( η g + η ( H + λ ) uitalic_t + η ϵitalic_t ) + μ ( uitalic_t - uitalic_t - 1 ) , t ≥ 1 , u- 1 := u0 Let H be symmetric and symmetric positive-definite. Suppose, μ∈(0,1)01μ∈(0,1)μ ∈ ( 0 , 1 ) and η<1/(λmax(H)+λ)1subscriptη<1/( _ (H)+λ)η < 1 / ( λroman_max ( H ) + λ ). Then, for u⋆=−(H+λ)−1gsuperscript⋆superscript1u =-(H+λ)^-1gu⋆ = - ( H + λ )- 1 g, ‖ut−u⋆‖≤(1−min1−μ,ηλ1−μ)t‖u0−u⋆‖+2maxη1−μ,1−μλmaxj≤t‖ϵj‖normsubscriptsuperscript⋆superscript111normsubscript0superscript⋆211subscriptnormsubscriptitalic-ϵ\|u_t-u \|≤ (1- \1- μ, ηλ1-% μ \ )^t\|u_0-u \|+ 2 \ η1-% μ, 1-μλ \ _j≤ t\| _j\|∥ uitalic_t - u⋆ ∥ ≤ ( 1 - min 1 - square-root start_ARG μ end_ARG , divide start_ARG η λ end_ARG start_ARG 1 - μ end_ARG )t ∥ u0 - u⋆ ∥ + square-root start_ARG 2 end_ARG max divide start_ARG η end_ARG start_ARG 1 - square-root start_ARG μ end_ARG end_ARG , divide start_ARG 1 - μ end_ARG start_ARG λ end_ARG maxitalic_j ≤ t ∥ ϵitalic_j ∥ Proof. Denote Hλ=H+λIsubscriptH_λ=H+λ IHitalic_λ = H + λ I. Let us rewrite the updates in the standard form, denoting the momentum vt=−(ut−ut−1)subscriptsubscriptsubscript1v_t=-(u_t-u_t-1)vitalic_t = - ( uitalic_t - uitalic_t - 1 ) we have v0=0subscript00v_0=0v0 = 0 and ut+1subscript1 u_t+1uitalic_t + 1 =(I−ηHλ)ut−μvt−ηg−ηϵtabsentsubscriptsubscriptsubscriptsubscriptitalic-ϵ =(I-η H_λ)u_t-μ v_t-η g-η _t= ( I - η Hitalic_λ ) uitalic_t - μ vitalic_t - η g - η ϵitalic_t vt+1subscript1 v_t+1vitalic_t + 1 =ηHλut+μvt+ηg+ηϵtabsentsubscriptsubscriptsubscriptsubscriptitalic-ϵ =η H_λu_t+μ v_t+η g+η _t= η Hitalic_λ uitalic_t + μ vitalic_t + η g + η ϵitalic_t Thus, we have the following linear update (ut+1vt+1)=(1−ηHλ−μηHλμ)(utvt)+(−ηgηg)+(−ηϵtηϵt)matrixsubscript1subscript1matrix1subscriptsubscriptmatrixsubscriptsubscriptmatrixmatrixsubscriptitalic-ϵsubscriptitalic-ϵ pmatrixu_t+1\\ v_t+1 pmatrix= pmatrix1-η H_λ&-μ\\ η H_λ&μ pmatrix pmatrixu_t\\ v_t pmatrix+ pmatrix-η g\\ η g pmatrix+ pmatrix-η _t\\ η _t pmatrix( start_ARG start_ROW start_CELL uitalic_t + 1 end_CELL end_ROW start_ROW start_CELL vitalic_t + 1 end_CELL end_ROW end_ARG ) = ( start_ARG start_ROW start_CELL 1 - η Hitalic_λ end_CELL start_CELL - μ end_CELL end_ROW start_ROW start_CELL η Hitalic_λ end_CELL start_CELL μ end_CELL end_ROW end_ARG ) ( start_ARG start_ROW start_CELL uitalic_t end_CELL end_ROW start_ROW start_CELL vitalic_t end_CELL end_ROW end_ARG ) + ( start_ARG start_ROW start_CELL - η g end_CELL end_ROW start_ROW start_CELL η g end_CELL end_ROW end_ARG ) + ( start_ARG start_ROW start_CELL - η ϵitalic_t end_CELL end_ROW start_ROW start_CELL η ϵitalic_t end_CELL end_ROW end_ARG ) (11) First, we want to show that the spectral norm of the matrix A:=(1−ηHλ−μηHλμ)assignmatrix1subscriptsubscriptA:= pmatrix1-η H_λ&-μ\\ η H_λ&μ pmatrixA := ( start_ARG start_ROW start_CELL 1 - η Hitalic_λ end_CELL start_CELL - μ end_CELL end_ROW start_ROW start_CELL η Hitalic_λ end_CELL start_CELL μ end_CELL end_ROW end_ARG ) is bounded by R:=maxμ,(1−ηλ1−μ)<1,assign111R:= \ μ, (1- ηλ1-μ ) \<1,R := max square-root start_ARG μ end_ARG , ( 1 - divide start_ARG η λ end_ARG start_ARG 1 - μ end_ARG ) < 1 , so that we can “unroll” the recursive equation. Bounding spectral norm. Given the decomposition Hλ=VDV⊤subscriptsuperscripttopH_λ=VDV Hitalic_λ = V D V⊤, we can rotate the matrix A turning it into into N independent updates with 2×2222× 22 × 2 matrices. Let D=diagh1,…,hNdiagsubscriptℎ1…subscriptℎD=diag\h_1,…,h_N\D = diag h1 , … , hitalic_N where h1≥⋯≥hN≥λsubscriptℎ1⋯subscriptℎh_1≥…≥ h_N≥ 1 ≥ ⋯ ≥ hitalic_N ≥ λ are the eigenvalues of H. By rotating the matrix A with (V00V)matrix00 pmatrixV&0\\ 0&V pmatrix( start_ARG start_ROW start_CELL V end_CELL start_CELL 0 end_CELL end_ROW start_ROW start_CELL 0 end_CELL start_CELL V end_CELL end_ROW end_ARG ), we can split into N non-intersecting 2×2222× 22 × 2 blocks on positions (i,N+i)×(i,N+i)(i,N+i)×(i,N+i)( i , N + i ) × ( i , N + i ). The corresponding 2×2222× 22 × 2 matrix updates look as follows (1−ηhi−μηhiμ)matrix1subscriptℎsubscriptℎ pmatrix1-η h_i&-μ\\ η h_i&μ pmatrix( start_ARG start_ROW start_CELL 1 - η hitalic_i end_CELL start_CELL - μ end_CELL end_ROW start_ROW start_CELL η hitalic_i end_CELL start_CELL μ end_CELL end_ROW end_ARG ) This matrix can have two real or two complex eigenvalues depending on the value hℎh (let us drop the index for now). Notice that the determinant equals (1−ηh)μ+ηhμ=μ1ℎ(1-η h)μ+η hμ=μ( 1 - η h ) μ + η h μ = μ, so in the latter case, the matrix is always a contraction, since both eigenvalues have norm μ μsquare-root start_ARG μ end_ARG. In the former case, we have two eigenvalues r1,2=1−ηh+μ2±(1−ηh+μ2)2−μsubscript12plus-or-minus1ℎ2superscript1ℎ22r_1,2= 1-η h+μ2± ( 1-η h+μ2 )^2% -μr1 , 2 = divide start_ARG 1 - η h + μ end_ARG start_ARG 2 end_ARG ± square-root start_ARG ( divide start_ARG 1 - η h + μ end_ARG start_ARG 2 end_ARG )2 - μ end_ARG Under the assumption ηh≤1ℎ1η h≤ 1η h ≤ 1, we have that both are positive and r2subscript2r_2r2 is the larger one. Since h≥λℎh≥ ≥ λ, we have that r2subscript2r_2r2 is at most r2≤1−ηλ+μ2+(1−ηλ+μ2)2−μ,subscript212superscript122r_2≤ 1-ηλ+μ2+ ( 1-ηλ+μ2% )^2-μ,r2 ≤ divide start_ARG 1 - η λ + μ end_ARG start_ARG 2 end_ARG + square-root start_ARG ( divide start_ARG 1 - η λ + μ end_ARG start_ARG 2 end_ARG )2 - μ end_ARG , and it is straightforward to check the following simplified bound r2≤1−ηλ1−μ.subscript211r_2≤ 1- ηλ1-μ.r2 ≤ 1 - divide start_ARG η λ end_ARG start_ARG 1 - μ end_ARG . Thus, the absolute maximum eigenvalue of the full matrix is bounded by R:=μ∨(1−ηλ1−μ)assign11R:= μ (1- ηλ1-μ )R := square-root start_ARG μ end_ARG ∨ ( 1 - divide start_ARG η λ end_ARG start_ARG 1 - μ end_ARG ). Expanding the recursion. Let us now expand the recursion in Eq. (11). Denoting A=(1−ηHλη(H+λ)μ)matrix1subscriptA= pmatrix1-η H_λ&η(H+λ)&μ pmatrixA = ( start_ARG start_ROW start_CELL 1 - η Hitalic_λ end_CELL start_CELL η ( H + λ ) end_CELL start_CELL μ end_CELL end_ROW end_ARG ), we have that (utvt)matrixsubscriptsubscript pmatrixu_t\\ v_t pmatrix( start_ARG start_ROW start_CELL uitalic_t end_CELL end_ROW start_ROW start_CELL vitalic_t end_CELL end_ROW end_ARG ) =A(ut−1vt−1)+(−ηgηg)+(−ηϵtηϵt)absentmatrixsubscript1subscript1matrixmatrixsubscriptitalic-ϵsubscriptitalic-ϵ =A pmatrixu_t-1\\ v_t-1 pmatrix+ pmatrix-η g\\ η g pmatrix+ pmatrix-η _t\\ η _t pmatrix= A ( start_ARG start_ROW start_CELL uitalic_t - 1 end_CELL end_ROW start_ROW start_CELL vitalic_t - 1 end_CELL end_ROW end_ARG ) + ( start_ARG start_ROW start_CELL - η g end_CELL end_ROW start_ROW start_CELL η g end_CELL end_ROW end_ARG ) + ( start_ARG start_ROW start_CELL - η ϵitalic_t end_CELL end_ROW start_ROW start_CELL η ϵitalic_t end_CELL end_ROW end_ARG ) =A2(ut−2vt−2)+(I+A)(−ηgηg)+(−ηϵtηϵt)+A(−ηϵt−1ηϵt−1)absentsuperscript2matrixsubscript2subscript2matrixmatrixsubscriptitalic-ϵsubscriptitalic-ϵmatrixsubscriptitalic-ϵ1subscriptitalic-ϵ1 =A^2 pmatrixu_t-2\\ v_t-2 pmatrix+(I+A) pmatrix-η g\\ η g pmatrix+ pmatrix-η _t\\ η _t pmatrix+A pmatrix-η _t-1\\ η _t-1 pmatrix= A2 ( start_ARG start_ROW start_CELL uitalic_t - 2 end_CELL end_ROW start_ROW start_CELL vitalic_t - 2 end_CELL end_ROW end_ARG ) + ( I + A ) ( start_ARG start_ROW start_CELL - η g end_CELL end_ROW start_ROW start_CELL η g end_CELL end_ROW end_ARG ) + ( start_ARG start_ROW start_CELL - η ϵitalic_t end_CELL end_ROW start_ROW start_CELL η ϵitalic_t end_CELL end_ROW end_ARG ) + A ( start_ARG start_ROW start_CELL - η ϵitalic_t - 1 end_CELL end_ROW start_ROW start_CELL η ϵitalic_t - 1 end_CELL end_ROW end_ARG ) =At(u00)+(I+⋯+At−1)(−ηgηg)+∑j=0t−1Aj(−ηϵt−jηϵt−j)absentsuperscriptmatrixsubscript00⋯superscript1matrixsuperscriptsubscript01superscriptmatrixsubscriptitalic-ϵsubscriptitalic-ϵ =A^t pmatrixu_0\\ 0 pmatrix+(I+…+A^t-1) pmatrix-η g\\ η g pmatrix+ _j=0^t-1A^j pmatrix-η _t-j\\ η _t-j pmatrix= Aitalic_t ( start_ARG start_ROW start_CELL u0 end_CELL end_ROW start_ROW start_CELL 0 end_CELL end_ROW end_ARG ) + ( I + ⋯ + Aitalic_t - 1 ) ( start_ARG start_ROW start_CELL - η g end_CELL end_ROW start_ROW start_CELL η g end_CELL end_ROW end_ARG ) + ∑j = 0t - 1 Aitalic_j ( start_ARG start_ROW start_CELL - η ϵitalic_t - j end_CELL end_ROW start_ROW start_CELL η ϵitalic_t - j end_CELL end_ROW end_ARG ) We bound the error term as ‖∑j=0t−1Aj(−ηϵt−jηϵt−j)‖≤∑j=0t−1Rj2η‖ϵt−j‖≤2η1−Rmaxj≤t‖ϵj‖.normsuperscriptsubscript01superscriptmatrixsubscriptitalic-ϵsubscriptitalic-ϵsuperscriptsubscript01superscript2normsubscriptitalic-ϵ21subscriptnormsubscriptitalic-ϵ \| _j=0^t-1A^j pmatrix-η _t-j\\ η _t-j pmatrix \|≤ _j=0^t-1R^j 2η\|% _t-j\|≤ 2η1-R _j≤ t\| _j\|\,.∥ ∑j = 0t - 1 Aitalic_j ( start_ARG start_ROW start_CELL - η ϵitalic_t - j end_CELL end_ROW start_ROW start_CELL η ϵitalic_t - j end_CELL end_ROW end_ARG ) ∥ ≤ ∑j = 0t - 1 Ritalic_j square-root start_ARG 2 end_ARG η ∥ ϵitalic_t - j ∥ ≤ divide start_ARG square-root start_ARG 2 end_ARG η end_ARG start_ARG 1 - R end_ARG maxitalic_j ≤ t ∥ ϵitalic_j ∥ . Let us now calculate (I−A)−1=(ηHλμηHλ1−μ)−1=((1−μ)(ηHλ)−1−μ(ηHλ)−1I)superscript1superscriptmatrixsubscriptsubscript11matrix1superscriptsubscript1superscriptsubscript1(I-A)^-1= pmatrixη H_λ&μ\\ η H_λ&1-μ pmatrix^-1= pmatrix(1-μ)(η H_% λ)^-1&-μ(η H_λ)^-1\\ I&I pmatrix( I - A )- 1 = ( start_ARG start_ROW start_CELL η Hitalic_λ end_CELL start_CELL μ end_CELL end_ROW start_ROW start_CELL η Hitalic_λ end_CELL start_CELL 1 - μ end_CELL end_ROW end_ARG )- 1 = ( start_ARG start_ROW start_CELL ( 1 - μ ) ( η Hitalic_λ )- 1 end_CELL start_CELL - μ ( η Hitalic_λ )- 1 end_CELL end_ROW start_ROW start_CELL I end_CELL start_CELL I end_CELL end_ROW end_ARG ) Then, we get that (I+⋯+At−1)(−ηgηg)⋯superscript1matrix (I+…+A^t-1) pmatrix-η g\\ η g pmatrix( I + ⋯ + Aitalic_t - 1 ) ( start_ARG start_ROW start_CELL - η g end_CELL end_ROW start_ROW start_CELL η g end_CELL end_ROW end_ARG ) =(I−At)(I−A)−1(−ηgηg)=(I−At)(−Hλ−1g0)absentsuperscriptsuperscript1matrixsuperscriptmatrixsuperscriptsubscript10 =(I-A^t)(I-A)^-1 pmatrix-η g\\ η g pmatrix=(I-A^t) pmatrix-H_λ^-1g\\ 0 pmatrix= ( I - Aitalic_t ) ( I - A )- 1 ( start_ARG start_ROW start_CELL - η g end_CELL end_ROW start_ROW start_CELL η g end_CELL end_ROW end_ARG ) = ( I - Aitalic_t ) ( start_ARG start_ROW start_CELL - Hitalic_λ- 1 g end_CELL end_ROW start_ROW start_CELL 0 end_CELL end_ROW end_ARG ) This brings us to the expression (ut+Hλ−1gvt)=At(u0+Hλ−1g0)+O≤(2η1−Rmaxj≤t‖ϵj‖)matrixsubscriptsuperscriptsubscript1subscriptsuperscriptmatrixsubscript0superscriptsubscript10subscript21subscriptnormsubscriptitalic-ϵ pmatrixu_t+H_λ^-1g\\ v_t pmatrix=A^t pmatrixu_0+H_λ^-1g\\ 0 pmatrix+O_≤ ( 2η1-R _j≤ t\| _% j\| )( start_ARG start_ROW start_CELL uitalic_t + Hitalic_λ- 1 g end_CELL end_ROW start_ROW start_CELL vitalic_t end_CELL end_ROW end_ARG ) = Aitalic_t ( start_ARG start_ROW start_CELL u0 + Hitalic_λ- 1 g end_CELL end_ROW start_ROW start_CELL 0 end_CELL end_ROW end_ARG ) + O≤ ( divide start_ARG square-root start_ARG 2 end_ARG η end_ARG start_ARG 1 - R end_ARG maxitalic_j ≤ t ∥ ϵitalic_j ∥ ) Given ‖At‖≤Rt=(1−min1−μ,ηλ1−μ)tnormsuperscriptsuperscriptsuperscript111\|A^t\|≤ R^t= (1- \1- μ, ηλ1-μ% \ )^t∥ Aitalic_t ∥ ≤ Ritalic_t = ( 1 - min 1 - square-root start_ARG μ end_ARG , divide start_ARG η λ end_ARG start_ARG 1 - μ end_ARG )t, we obtain the statement of the lemma. ∎ D.2 Remaining proof of Theorem 3.1 Notation. v=O≤(x)subscriptv=O_≤(x)v = O≤ ( x ) means ‖v‖≤xnorm\|v\|≤ x∥ v ∥ ≤ x. This notation is more convenient for proving inequalities by induction. We also denote Hλ(θ)=H(θ)+λIsubscriptH_λ(θ)=H(θ)+λ IHitalic_λ ( θ ) = H ( θ ) + λ I for all λ, and Nλ(θ)=Hλ(θ)−1∇L(θ)subscriptsubscriptsuperscript1∇N_λ(θ)=H_λ(θ)^-1∇ L(θ)Nitalic_λ ( θ ) = Hitalic_λ ( θ )- 1 ∇ L ( θ ) is the natural gradient. We prove the theorem under the following regularity conditions. Condition D.2. There is a constant K≥11K≥ 1K ≥ 1 such that for any θ,θ′θ,θ θ , θ′ and λ′∈λ,λ¯superscript′¯λ ∈\λ, λ\λ′ ∈ λ , over¯ start_ARG λ end_ARG , ‖Hλ′(θ)‖normsubscriptsuperscript′ \|H_λ (θ)\|∥ Hitalic_λ′ ( θ ) ∥ ≤K,absent ≤ K,≤ K , ‖Hλ′−1(θ)∇L(θ)−Hλ′(θ′)−1∇L(θ)‖normsuperscriptsubscriptsuperscript′1∇subscriptsuperscript′superscript′1∇ \|H_λ ^-1(θ)∇ L(θ)-H_λ^% (θ )^-1∇ L(θ)\|∥ Hitalic_λ′- 1 ( θ ) ∇ L ( θ ) - Hitalic_λ′ ( θ′ )- 1 ∇ L ( θ ) ∥ ≤K‖θ−θ′‖,absentnormsuperscript′ ≤ K\|θ-θ \|,≤ K ∥ θ - θ′ ∥ , ‖∇θ(θ,θ′)−H(θ)(θ−θ′)‖normsubscript∇superscript′ \| _θD(θ,θ )-H(θ)(% θ-θ )\|∥ ∇θ D ( θ , θ′ ) - H ( θ ) ( θ - θ′ ) ∥ ≤K‖θ−θ′‖2,absentsuperscriptnormsuperscript′2 ≤ K\|θ-θ \|^2,≤ K ∥ θ - θ′ ∥2 , ‖Hλ′−1(θ)∇L(θ)‖normsuperscriptsubscriptsuperscript′1∇ \|H_λ ^-1(θ)∇ L(θ)\|∥ Hitalic_λ′- 1 ( θ ) ∇ L ( θ ) ∥ ≤K.absent ≤ K\,.≤ K . The third condition quantifies, how well the local quadratic approximation in Eq. 5 holds. Define ut=θt−θt′subscriptsubscriptsuperscriptsubscript′u_t= _t- _t uitalic_t = θitalic_t - θitalic_t′. Furthermore, define Ut:=maxt′≤t‖ut′‖,Dt:=maxt′≤t‖ut′−ut′−1‖.formulae-sequenceassignsubscriptsubscriptsuperscript′normsubscriptsuperscript′assignsubscriptsubscriptsuperscript′normsubscriptsuperscript′subscriptsuperscript′1U_t:= _t ≤ t\|u_t \|, D_t:= _t % ≤ t\|u_t -u_t -1\|\,.Uitalic_t := maxitalic_t′ ≤ t ∥ uitalic_t′ ∥ , Ditalic_t := maxitalic_t′ ≤ t ∥ uitalic_t′ - uitalic_t′ - 1 ∥ . (12) Let us rewrite the recursion θt+1−θt−μ(θt−θt−1)subscript1subscriptsubscriptsubscript1 _t+1- _t-μ( _t- _t-1)θitalic_t + 1 - θitalic_t - μ ( θitalic_t - θitalic_t - 1 ) =−η[α∇L(θt)+∇θ(θt,θt′)+λ(θt−θt′)]absentdelimited-[]∇subscriptsubscript∇subscriptsuperscriptsubscript′subscriptsuperscriptsubscript′ =-η [α∇ L( _t)+ _θD% ( _t, _t )+λ( _t- _t ) ]= - η [ α ∇ L ( θitalic_t ) + ∇θ D ( θitalic_t , θitalic_t′ ) + λ ( θitalic_t - θitalic_t′ ) ] =−η[α∇L(θt)+H(θt)(θt−θt′)+λ(θt−θt′)]+KηO≤(‖θt−θt′‖2)absentdelimited-[]∇subscriptsubscriptsubscriptsuperscriptsubscript′subscriptsuperscriptsubscript′subscriptsuperscriptnormsubscriptsuperscriptsubscript′2 =-η [α∇ L( _t)+H( _t)( _t-% _t )+λ( _t- _t ) ]+Kη\,O_% ≤(\| _t- _t \|^2)= - η [ α ∇ L ( θitalic_t ) + H ( θitalic_t ) ( θitalic_t - θitalic_t′ ) + λ ( θitalic_t - θitalic_t′ ) ] + K η O≤ ( ∥ θitalic_t - θitalic_t′ ∥2 ) =−η[α∇L(θt)+H(θt)(θt−θt′)+λ(θt−θt′)]+KηO≤(Ut2).absentdelimited-[]∇subscriptsubscriptsubscriptsuperscriptsubscript′subscriptsuperscriptsubscript′subscriptsuperscriptsubscript2 =-η [α∇ L( _t)+H( _t)( _t-% _t )+λ( _t- _t ) ]+Kη\,O_% ≤(U_t^2)\,.= - η [ α ∇ L ( θitalic_t ) + H ( θitalic_t ) ( θitalic_t - θitalic_t′ ) + λ ( θitalic_t - θitalic_t′ ) ] + K η O≤ ( Uitalic_t2 ) . Therefore, θt+1−θt−μ(θt−θt−1)=−η[α∇L(θt)+Hλ(θt)ut]+KηO≤(Ut2).subscript1subscriptsubscriptsubscript1delimited-[]∇subscriptsubscriptsubscriptsubscriptsubscriptsuperscriptsubscript2 _t+1- _t-μ( _t- _t-1)=-η [α∇ L% ( _t)+H_λ( _t)u_t ]+Kη O_≤(U_t^2).θitalic_t + 1 - θitalic_t - μ ( θitalic_t - θitalic_t - 1 ) = - η [ α ∇ L ( θitalic_t ) + Hitalic_λ ( θitalic_t ) uitalic_t ] + K η O≤ ( Uitalic_t2 ) . (13) We want to construct a similar update for sequence utsubscriptu_tuitalic_t, for which we need to subtract the corresponding combination θt+1′−θt′−μ(θt′−θt−1′)superscriptsubscript1′subscript′subscript′subscript1′ _t+1 - _t -μ( _t - _t-1% )θitalic_t + 1′ - θitalic_t′ - μ ( θitalic_t′ - θitalic_t - 1′ ). Using the sliding reference update and the substitution θt+1=θt+1′+ut+1subscript1superscriptsubscript1′subscript1 _t+1= _t+1 +u_t+1θitalic_t + 1 = θitalic_t + 1′ + uitalic_t + 1, we have that θt+1′subscript1′ _t+1 θitalic_t + 1′ =(1−ηκ)θt′+ηκ(θt+1′+ut+1),absent1superscriptsubscript′subscript1′subscript1 =(1-ηκ) _t +ηκ( _t+1^% +u_t+1),= ( 1 - η κ ) θitalic_t′ + η κ ( θitalic_t + 1′ + uitalic_t + 1 ) , (1−ηκ)(θt+1′−θt′)1superscriptsubscript1′subscript′ (1-ηκ)( _t+1 - _t )( 1 - η κ ) ( θitalic_t + 1′ - θitalic_t′ ) =ηκut+1,absentsubscript1 =ηκ u_t+1,= η κ uitalic_t + 1 , θt+1′−θt′subscript1′subscript′ _t+1 - _t θitalic_t + 1′ - θitalic_t′ =ηκ¯ut+1,absent¯subscript1 =η κu_t+1,= η over¯ start_ARG κ end_ARG uitalic_t + 1 , θt+1′−θt′−μ(θt′−θt−1′)superscriptsubscript1′subscript′subscript′subscript1′ _t+1 - _t -μ( _t % - _t-1 )θitalic_t + 1′ - θitalic_t′ - μ ( θitalic_t′ - θitalic_t - 1′ ) =ηκ¯(ut+1−ut)+(1−μ)ηκ¯ut,absent¯subscript1subscript1¯subscript =η κ(u_t+1-u_t)+(1-μ)η % κu_t,= η over¯ start_ARG κ end_ARG ( uitalic_t + 1 - uitalic_t ) + ( 1 - μ ) η over¯ start_ARG κ end_ARG uitalic_t , where we define κ¯=κ/(1−ηκ)¯1 κ=κ/(1-ηκ)over¯ start_ARG κ end_ARG = κ / ( 1 - η κ ). Subtracting the latter equation from (13), we obtain an iterative update for sequence utsubscriptu_tuitalic_t, (1+ηκ¯)(ut+1−ut)−μ(ut−ut−1)=−η[α∇L(θt)+(Hλ(θt)+(1−μ)κ¯I)ut]+KηO≤(Ut2).1¯subscript1subscriptsubscriptsubscript1delimited-[]∇subscriptsubscriptsubscript1¯subscriptsubscriptsuperscriptsubscript2(1+η κ)(u_t+1-u_t)-μ(u_t-u_t-1)=-η [α% ∇ L( _t)+(H_λ( _t)+(1-μ) κI)u_t% ]+Kη O_≤(U_t^2)\,.( 1 + η over¯ start_ARG κ end_ARG ) ( uitalic_t + 1 - uitalic_t ) - μ ( uitalic_t - uitalic_t - 1 ) = - η [ α ∇ L ( θitalic_t ) + ( Hitalic_λ ( θitalic_t ) + ( 1 - μ ) over¯ start_ARG κ end_ARG I ) uitalic_t ] + K η O≤ ( Uitalic_t2 ) . We next divide both sides by 1+ηκ¯=1/(1−ηκ)1¯111+η κ=1/(1-ηκ)1 + η over¯ start_ARG κ end_ARG = 1 / ( 1 - η κ ), and simplify the notation with λ¯=λ+(1−μ)κ¯1¯ λ=λ+(1-μ) κover¯ start_ARG λ end_ARG = λ + ( 1 - μ ) over¯ start_ARG κ end_ARG, η¯=η(1−ηκ)¯1 η=η(1-ηκ)over¯ start_ARG η end_ARG = η ( 1 - η κ ), and μ¯=μ(1−ηκ)¯1 μ=μ(1-ηκ)over¯ start_ARG μ end_ARG = μ ( 1 - η κ ), which leads us to the following recursive updates, ut+1−ut−μ¯(ut−ut−1)=−η¯[α∇L(θt)+Hλ¯(θt)ut]+KηO≤(Ut2).subscript1subscript¯subscriptsubscript1¯delimited-[]∇subscriptsubscript¯subscriptsubscriptsubscriptsuperscriptsubscript2u_t+1-u_t- μ(u_t-u_t-1)=- η [α∇ L% ( _t)+H_ λ( _t)u_t ]+Kη O_≤(U_% t^2)\,.uitalic_t + 1 - uitalic_t - over¯ start_ARG μ end_ARG ( uitalic_t - uitalic_t - 1 ) = - over¯ start_ARG η end_ARG [ α ∇ L ( θitalic_t ) + Hover¯ start_ARG λ end_ARG ( θitalic_t ) uitalic_t ] + K η O≤ ( Uitalic_t2 ) . (14) Notice that these updates look similar to the IHVP updates in Lemma D.1. Here, θtsubscript _tθitalic_t slightly changes with each update as well, and correspondingly the gradient ∇L(θt)∇subscript∇ L( _t)∇ L ( θitalic_t ) and the Hessian Hλ¯(θt)subscript¯subscriptH_ λ( _t)Hover¯ start_ARG λ end_ARG ( θitalic_t ). We account for it through the error term ϵtsubscriptitalic-ϵ _tϵitalic_t in the lemma, and for sufficiently small α, this is not a significant change in parameter and eventually after certain amount of steps we reach approximation ut≈−αHλ¯(θt)−1∇L(θt)subscriptsubscript¯superscriptsubscript1∇subscriptu_t≈-α H_ λ( _t)^-1∇ L( _t)uitalic_t ≈ - α Hover¯ start_ARG λ end_ARG ( θitalic_t )- 1 ∇ L ( θitalic_t ) along the remaining training trajectory. This will allow us to show that the trajectory θtsubscript _tθitalic_t matches the natural gradient trajectory. Indeed, using the sliding reference equation and the substitution θt′=θt−utsuperscriptsubscript′subscriptsubscript _t = _t-u_tθitalic_t′ = θitalic_t - uitalic_t, we have θt+1−ut+1=(1−κη)(θt−ut)+κηθt+1,subscript1subscript11subscriptsubscriptsubscript1 _t+1-u_t+1=(1-κη)( _t-u_t)+κη _t+1,θitalic_t + 1 - uitalic_t + 1 = ( 1 - κ η ) ( θitalic_t - uitalic_t ) + κ η θitalic_t + 1 , which after some rearrangements turns into θt+1=θt+κ¯ηut+1+(ut+1−ut).subscript1subscript¯subscript1subscript1subscript _t+1= _t+ κη u_t+1+(u_t+1-u_t)\,.θitalic_t + 1 = θitalic_t + over¯ start_ARG κ end_ARG η uitalic_t + 1 + ( uitalic_t + 1 - uitalic_t ) . (15) Given that utsubscriptu_tuitalic_t approximates the natural gradient −αHλ+κ¯(θt)−1∇L(θt)subscript¯superscriptsubscript1∇subscript-α H_λ+ κ( _t)^-1∇ L( _t)- α Hitalic_λ + over¯ start_ARG κ end_ARG ( θitalic_t )- 1 ∇ L ( θitalic_t ), we can see that this update looks more like natural gradient descent, granted that the difference ut−ut+1subscriptsubscript1u_t-u_t+1uitalic_t - uitalic_t + 1 has smaller order. To complete the proof we are going to do the following steps: 1. Show that Ut=O(α)subscriptU_t=O(α)Uitalic_t = O ( α ) and Dt=O(ηα2)subscriptsuperscript2D_t=O(ηα^2)Ditalic_t = O ( η α2 ) for t≤O(1/(ηα))1t≤ O(1/(ηα))t ≤ O ( 1 / ( η α ) ) by induction. 2. Show that ut≈−αNλ¯(θt)subscriptsubscript¯subscriptu_t≈-α N_ λ( _t)uitalic_t ≈ - α Nover¯ start_ARG λ end_ARG ( θitalic_t ) in an early snapshot t∈[t0,2t0]subscript02subscript0t∈[t_0,2t_0]t ∈ [ t0 , 2 t0 ], where t0=O(1ηlog(1/α))subscript011t_0=O( 1η (1/α))t0 = O ( divide start_ARG 1 end_ARG start_ARG η end_ARG log ( 1 / α ) ). 3. Show that this approximation is maintained for the remaining steps, simultaneously maintaining the bounds on UtsubscriptU_tUitalic_t, DtsubscriptD_tDitalic_t. This is done by induction as well, but with a different induction assumption. 4. Connect θtsubscript _tθitalic_t to a natural gradient trajectory using Eq. (15). Stage 1: convergence to NG. At the start of iterations, the sequence utsubscriptu_tuitalic_t is initialized at 00, since at first we have that θ0=θ0′subscript0superscriptsubscript0′ _0= _0 θ0 = θ0′, and it does not approximate the natural gradient. According to Lemma D.1, the IHVP updates converge in approximately t0=Ω((1−μ)/(λη))subscript0Ω1t_0= ((1-μ)/(λη))t0 = Ω ( ( 1 - μ ) / ( λ η ) ) steps, and our goal is to show that about the same number of steps is sufficient for the iterations (14). Because we do not expect the parameter θ to change dramatically within this number of steps, let us compare everything at parameter θ0subscript0 _0θ0. From (15), we have that ‖θt−θt−1‖≤Atnormsubscriptsubscript1subscript\| _t- _t-1\|≤ A_t∥ θitalic_t - θitalic_t - 1 ∥ ≤ Aitalic_t, where we set At:=κ¯ηUt+Dtassignsubscript¯subscriptsubscriptA_t:= κη U_t+D_tAitalic_t := over¯ start_ARG κ end_ARG η Uitalic_t + Ditalic_t. Therefore, we also have ‖θt−θ0‖≤tAtnormsubscriptsubscript0subscript\| _t- _0\|≤ tA_t∥ θitalic_t - θ0 ∥ ≤ t Aitalic_t. Then, our regularity conditions yield ∇L(θt)=∇L(θ0)+O≤(KtAt),Hλ¯(θt)ut=Hλ¯(θ0)ut+O≤(KtAtUt).formulae-sequence∇subscript∇subscript0subscriptsubscriptsubscript¯subscriptsubscriptsubscript¯subscript0subscriptsubscriptsubscriptsubscript∇ L( _t)=∇ L( _0)+O_≤ (KtA_t ), H% _ λ( _t)u_t=H_ λ( _0)u_t% +O_≤(KtA_tU_t).∇ L ( θitalic_t ) = ∇ L ( θ0 ) + O≤ ( K t Aitalic_t ) , Hover¯ start_ARG λ end_ARG ( θitalic_t ) uitalic_t = Hover¯ start_ARG λ end_ARG ( θ0 ) uitalic_t + O≤ ( K t Aitalic_t Uitalic_t ) . The equation (14) can therefore be approximated as follows ut+1−ut−μ¯(ut−ut−1)=−η¯[α∇L(θ0)+Hλ¯(θ0)ut]+η¯O≤(KtAt(α+Ut)+Ut2).subscript1subscript¯subscriptsubscript1¯delimited-[]∇subscript0subscript¯subscript0subscript¯subscriptsubscriptsubscriptsuperscriptsubscript2u_t+1-u_t- μ(u_t-u_t-1)=- η [α∇ L% ( _0)+H_ λ( _0)u_t ]+ ηO_% ≤ (KtA_t(α+U_t)+U_t^2 ).uitalic_t + 1 - uitalic_t - over¯ start_ARG μ end_ARG ( uitalic_t - uitalic_t - 1 ) = - over¯ start_ARG η end_ARG [ α ∇ L ( θ0 ) + Hover¯ start_ARG λ end_ARG ( θ0 ) uitalic_t ] + over¯ start_ARG η end_ARG O≤ ( K t Aitalic_t ( α + Uitalic_t ) + Uitalic_t2 ) . (16) We apply Lemma D.1 with ϵt=η¯O≤(KtAt(α+Ut)+Ut2)subscriptitalic-ϵ¯subscriptsubscriptsubscriptsuperscriptsubscript2 _t= ηO_≤ (KtA_t(α+U_t)+U_t^2 )ϵitalic_t = over¯ start_ARG η end_ARG O≤ ( K t Aitalic_t ( α + Uitalic_t ) + Uitalic_t2 ). Notice that in the above equation the momentum coefficient is replaced by μ¯=μ(1+ηκ)¯1 μ=μ(1+ηκ)over¯ start_ARG μ end_ARG = μ ( 1 + η κ ), which we assume is less than 1111. Denote for convenience, δ:=min1−μ¯,η¯λ¯(1−μ¯)=η¯λ¯(1−μ¯)assign1¯1¯1¯δ:= \1- μ, η % λ(1- μ) \= η λ% (1- μ)δ := min 1 - square-root start_ARG over¯ start_ARG μ end_ARG end_ARG , divide start_ARG over¯ start_ARG η end_ARG over¯ start_ARG λ end_ARG end_ARG start_ARG ( 1 - over¯ start_ARG μ end_ARG ) end_ARG = divide start_ARG over¯ start_ARG η end_ARG over¯ start_ARG λ end_ARG end_ARG start_ARG ( 1 - over¯ start_ARG μ end_ARG ) end_ARG where we us assume that η is sufficiently small, and we also assume that maxη¯1−μ¯,1−μ¯λ¯=(1−μ¯)/λ¯1¯1¯1¯ \ η1- μ, 1- % μ λ \=(1- μ)/ λmax divide start_ARG over¯ start_ARG η end_ARG end_ARG start_ARG 1 - square-root start_ARG over¯ start_ARG μ end_ARG end_ARG end_ARG , divide start_ARG 1 - over¯ start_ARG μ end_ARG end_ARG start_ARG over¯ start_ARG λ end_ARG end_ARG = ( 1 - over¯ start_ARG μ end_ARG ) / over¯ start_ARG λ end_ARG for the same reason. We first bound Ut,DtsubscriptsubscriptU_t,D_tUitalic_t , Ditalic_t in the initial steps by applying Lemma D.1 with a trivial inequality (1−δ)t<1superscript11(1-δ)^t<1( 1 - δ )t < 1, by induction. Lemma D.3. There is a c0=c0(K,κ,λ,μ)subscript0subscript0c_0=c_0(K,κ,λ,μ)c0 = c0 ( K , κ , λ , μ ) such that for t≤c0/(αη)subscript0t≤ c_0/(αη)t ≤ c0 / ( α η ) and sufficiently small η, Ut≤3Kα,Dt≤5(1−μ¯)−1K2ηα.formulae-sequencesubscript3subscript5superscript1¯1superscript2U_t≤ 3Kα, D_t≤ 5(1- μ)^-1K^2ηα.Uitalic_t ≤ 3 K α , Ditalic_t ≤ 5 ( 1 - over¯ start_ARG μ end_ARG )- 1 K2 η α . Proof. First we use a trivial inequality (1−δ)t<1superscript11(1-δ)^t<1( 1 - δ )t < 1 and the fact that ‖Hλ¯(θ)−1∇L(θ)‖≤Knormsubscript¯superscript1∇\|H_ λ(θ)^-1∇ L(θ)\|≤ K∥ Hover¯ start_ARG λ end_ARG ( θ )- 1 ∇ L ( θ ) ∥ ≤ K to obtain, ‖ut+1+αNλ¯(θ0)‖<α∥Nλ¯(θ0)∥+21−μ¯λ¯(KtAt(α+Ut)+Ut2)conditionalnormsubscript1subscript¯subscript0brasubscript¯subscript021¯subscriptsubscriptsuperscriptsubscript2\|u_t+1+α N_ λ( _0)\|<α\|N_ % λ( _0)\|+ 2 1- μ λ% (KtA_t(α+U_t)+U_t^2 )∥ uitalic_t + 1 + α Nover¯ start_ARG λ end_ARG ( θ0 ) ∥ < α ∥ Nover¯ start_ARG λ end_ARG ( θ0 ) ∥ + square-root start_ARG 2 end_ARG divide start_ARG 1 - over¯ start_ARG μ end_ARG end_ARG start_ARG over¯ start_ARG λ end_ARG end_ARG ( K t Aitalic_t ( α + Uitalic_t ) + Uitalic_t2 ) This first gives us ‖ut+1‖≤2Kα+2(1−μ¯)λ¯−1(KtAt(α+Ut)+Ut2).normsubscript1221¯superscript¯1subscriptsubscriptsuperscriptsubscript2\|u_t+1\|≤ 2Kα+ 2(1- μ) λ^-1(% KtA_t(α+U_t)+U_t^2).∥ uitalic_t + 1 ∥ ≤ 2 K α + square-root start_ARG 2 end_ARG ( 1 - over¯ start_ARG μ end_ARG ) over¯ start_ARG λ end_ARG- 1 ( K t Aitalic_t ( α + Uitalic_t ) + Uitalic_t2 ) . We then have by induction assumptions Ut≤3Kαsubscript3U_t≤ 3K _t ≤ 3 K α, Dt≤5(1−μ¯)−1K2ηαsubscript5superscript1¯1superscript2D_t≤ 5(1- μ)^-1K^2η _t ≤ 5 ( 1 - over¯ start_ARG μ end_ARG )- 1 K2 η α Ut+1subscript1 U_t+1Uitalic_t + 1 ≤2Kα+2(1−μ¯)λ¯−1[Kt(κ¯ηUt+Dt)(α+Ut)+Ut2]absent221¯superscript¯1delimited-[]¯subscriptsubscriptsubscriptsuperscriptsubscript2 ≤ 2Kα+ 2(1- μ) λ^-1% [Kt( κη U_t+D_t)(α+U_t)+U_t^2]≤ 2 K α + square-root start_ARG 2 end_ARG ( 1 - over¯ start_ARG μ end_ARG ) over¯ start_ARG λ end_ARG- 1 [ K t ( over¯ start_ARG κ end_ARG η Uitalic_t + Ditalic_t ) ( α + Uitalic_t ) + Uitalic_t2 ] ≤2Kα+2(1−μ¯)λ¯−1[Kt(κ¯η(3Kα)+5(1−μ)−1K2ηα)(3K+1)α+(3Kα)2]absent221¯superscript¯1delimited-[]¯35superscript11superscript231superscript32 ≤ 2Kα+ 2(1- μ) λ^-1% [Kt( κη(3Kα)+5(1-μ)^-1K^2ηα)(3K+1)α+% (3Kα)^2]≤ 2 K α + square-root start_ARG 2 end_ARG ( 1 - over¯ start_ARG μ end_ARG ) over¯ start_ARG λ end_ARG- 1 [ K t ( over¯ start_ARG κ end_ARG η ( 3 K α ) + 5 ( 1 - μ )- 1 K2 η α ) ( 3 K + 1 ) α + ( 3 K α )2 ] =2Kα+Kα×2(1−μ¯)(λ¯)−1(tηα3κ¯K+5(1−μ)−1K2(3K+1)+9Kα),absent221¯superscript¯13¯5superscript11superscript2319 =2Kα+Kα× 2(1- μ)( % λ)^-1 (tηα\3 κK+5(1-μ)^-1K^2\(3K+1% )+9Kα ),= 2 K α + K α × square-root start_ARG 2 end_ARG ( 1 - over¯ start_ARG μ end_ARG ) ( over¯ start_ARG λ end_ARG )- 1 ( t η α 3 over¯ start_ARG κ end_ARG K + 5 ( 1 - μ )- 1 K2 ( 3 K + 1 ) + 9 K α ) , with the last expression at most 3Kα33K 3 K α in the case where α≤λ¯/(182K(1−μ¯))¯1821¯α≤ λ/(18 2K(1- μ))α ≤ over¯ start_ARG λ end_ARG / ( 18 square-root start_ARG 2 end_ARG K ( 1 - over¯ start_ARG μ end_ARG ) ) and t≤c0/(ηα)subscript0t≤ c_0/(ηα)t ≤ c0 / ( η α ) with c0=λ¯(223(1−μ¯)κ¯K+5K2(3K+1))−1subscript0¯superscript2231¯5superscript2311c_0= λ (2 2\3(1- μ) κK% +5K^2\(3K+1) )^-1c0 = over¯ start_ARG λ end_ARG ( 2 square-root start_ARG 2 end_ARG 3 ( 1 - over¯ start_ARG μ end_ARG ) over¯ start_ARG κ end_ARG K + 5 K2 ( 3 K + 1 ) )- 1. We allow it to depend on κ¯ κover¯ start_ARG κ end_ARG since for small enough η it is at most 2κ22 2 κ. And then we also get from Eq. (16), ‖ut+1−ut‖≤μ¯‖ut−ut−1‖+η¯K(α+Ut)+η¯(KtAt(α+Ut)+Ut2).normsubscript1subscript¯normsubscriptsubscript1¯subscript¯subscriptsubscriptsuperscriptsubscript2\|u_t+1-u_t\|≤ μ\|u_t-u_t-1\|+ ηK(α+U% _t)+ η(KtA_t(α+U_t)+U_t^2).∥ uitalic_t + 1 - uitalic_t ∥ ≤ over¯ start_ARG μ end_ARG ∥ uitalic_t - uitalic_t - 1 ∥ + over¯ start_ARG η end_ARG K ( α + Uitalic_t ) + over¯ start_ARG η end_ARG ( K t Aitalic_t ( α + Uitalic_t ) + Uitalic_t2 ) . Then, Dt+1subscript1 D_t+1Ditalic_t + 1 ≤μ¯Dt+ηK(α+Ut)+η(KtAt(α+Ut)+Ut2)absent¯subscriptsubscriptsubscriptsubscriptsuperscriptsubscript2 ≤ μD_t+η K(α+U_t)+η(KtA_t(α+% U_t)+U_t^2)≤ over¯ start_ARG μ end_ARG Ditalic_t + η K ( α + Uitalic_t ) + η ( K t Aitalic_t ( α + Uitalic_t ) + Uitalic_t2 ) ≤5μ¯(1−μ¯)−1K2ηα+ηαK(3K+1)+η(Kt(κ¯η3Kα+5(1−μ¯)−1K2ηα)(3K+1)α+9K2α2)absent5¯superscript1¯1superscript231¯35superscript1¯1superscript2319superscript2superscript2 ≤ 5 μ(1- μ)^-1K^2ηα+η% α K(3K+1)+η (Kt( κη 3Kα+5(1- μ)% ^-1K^2ηα)(3K+1)α+9K^2α^2 )≤ 5 over¯ start_ARG μ end_ARG ( 1 - over¯ start_ARG μ end_ARG )- 1 K2 η α + η α K ( 3 K + 1 ) + η ( K t ( over¯ start_ARG κ end_ARG η 3 K α + 5 ( 1 - over¯ start_ARG μ end_ARG )- 1 K2 η α ) ( 3 K + 1 ) α + 9 K2 α2 ) =(1−μ¯)−1K2ηα5μ¯+(3+1/K)(1−μ¯)+(3(1−μ¯)+5K)(3K+1)αηt+9α,absentsuperscript1¯1superscript25¯311¯31¯5319 =(1- μ)^-1K^2ηα \5 μ+(3+% 1/K)(1- μ)+(3(1- μ)+5K)(3K+1)αη t+9α% \,= ( 1 - over¯ start_ARG μ end_ARG )- 1 K2 η α 5 over¯ start_ARG μ end_ARG + ( 3 + 1 / K ) ( 1 - over¯ start_ARG μ end_ARG ) + ( 3 ( 1 - over¯ start_ARG μ end_ARG ) + 5 K ) ( 3 K + 1 ) α η t + 9 α , with the value in the brackets at most 5555 as long as α≤(1−μ¯)/181¯18α≤(1- μ)/18α ≤ ( 1 - over¯ start_ARG μ end_ARG ) / 18 and t≤c0/(αη)subscript0t≤ c_0/(αη)t ≤ c0 / ( α η ) for c0=(1−μ¯)((6(1−μ¯)+10K)(3K+1))−1subscript01¯superscript61¯10311c_0=(1- μ) ((6(1- μ)+10K)(3K+1) )^-1c0 = ( 1 - over¯ start_ARG μ end_ARG ) ( ( 6 ( 1 - over¯ start_ARG μ end_ARG ) + 10 K ) ( 3 K + 1 ) )- 1, where we also use that K≥11K≥ 1K ≥ 1. ∎ Notice that t=c0/(αη)subscript0t=c_0/(αη)t = c0 / ( α η ) should be enough to have a reasonable approximation of the natural gradient, since it can be much bigger than (1−μ)/(λη)1(1-μ)/(λη)( 1 - μ ) / ( λ η ) when α is sufficiently small. Next, we use Lemma D.1 again. Now that we have At≤(3Kκ¯+5(1−μ¯)−1K2)αηsubscript3¯5superscript1¯1superscript2A_t≤(3K κ+5(1- μ)^-1K^2)α _t ≤ ( 3 K over¯ start_ARG κ end_ARG + 5 ( 1 - over¯ start_ARG μ end_ARG )- 1 K2 ) α η, we get that ∥ut+αNλ¯(θ0)∥≤K(1−η¯λ¯/(1−μ¯))t+λ¯−1K2((3κ¯+5K(1−μ¯)−1(tη)α2+9α2)\|u_t+α N_ λ( _0)\|≤ K(1- η% λ/(1- μ))^t+ λ^-1K^2 ((% 3 κ+5K(1- μ)^-1(tη)α^2+9α^2 )∥ uitalic_t + α Nover¯ start_ARG λ end_ARG ( θ0 ) ∥ ≤ K ( 1 - over¯ start_ARG η end_ARG over¯ start_ARG λ end_ARG / ( 1 - over¯ start_ARG μ end_ARG ) )t + over¯ start_ARG λ end_ARG- 1 K2 ( ( 3 over¯ start_ARG κ end_ARG + 5 K ( 1 - over¯ start_ARG μ end_ARG )- 1 ( t η ) α2 + 9 α2 ) Let us take t that satisfies (1−ηλ)t=α2superscript1superscript2(1-ηλ)^t=α^2( 1 - η λ )t = α2, so that the two terms have approximately the same order, which is t0=2(1−μ¯)η¯λ¯⌈log1α⌉subscript021¯1t_0= 2(1- μ) η λ % 1α 0 = divide start_ARG 2 ( 1 - over¯ start_ARG μ end_ARG ) end_ARG start_ARG over¯ start_ARG η end_ARG over¯ start_ARG λ end_ARG end_ARG ⌈ log divide start_ARG 1 end_ARG start_ARG α end_ARG ⌉. For sufficiently small α, the whole range [t0,2t0]subscript02subscript0[t_0,2t_0][ t0 , 2 t0 ] lies within c0/(ηα)subscript0c_0/(ηα)c0 / ( η α ), so that the bounds from Lemma D.3 hold. This yields for some C0=C0(K,κ,λ,μ)subscript0subscript0C_0=C_0(K,κ,λ,μ)C0 = C0 ( K , κ , λ , μ ) and sufficiently small α,ηα,ηα , η maxt∈[t0,2t0]‖ut+αHλ+κ¯(θt)−1∇L(θt)‖≤C0α2log1α.subscriptsubscript02subscript0normsubscriptsubscript¯superscriptsubscript1∇subscriptsubscript0superscript21 _t∈[t_0,2t_0]\|u_t+α H_λ+ κ( _% t)^-1∇ L( _t)\|≤ C_0α^2 1α\,.maxitalic_t ∈ [ t start_POSTSUBSCRIPT 0 , 2 t0 ] end_POSTSUBSCRIPT ∥ uitalic_t + α Hitalic_λ + over¯ start_ARG κ end_ARG ( θitalic_t )- 1 ∇ L ( θitalic_t ) ∥ ≤ C0 α2 log divide start_ARG 1 end_ARG start_ARG α end_ARG . We will assume that the choice of C0subscript0C_0C0 also ensures this inequality for t0+1subscript01t_0+1t0 + 1 in place of t0subscript0t_0t0. Stage 2: maintaining NG approximation. Next, we show that we can maintain the approximation of the natural gradient by utsubscriptu_tuitalic_t after step t≥t0subscript0t≥ t_0t ≥ t0 as well. We are going to use a stronger inequality to carry through by induction. We first slightly refine the bound on ut−ut−1subscriptsubscript1u_t-u_t-1uitalic_t - uitalic_t - 1 for t=t0,t0+1subscript0subscript01t=t_0,t_0+1t = t0 , t0 + 1. Since ‖θt−θt−1‖≤At≤C1ηαnormsubscriptsubscript1subscriptsubscript1\| _t- _t-1\|≤ A_t≤ C_1ηα∥ θitalic_t - θitalic_t - 1 ∥ ≤ Aitalic_t ≤ C1 η α, we get ‖ut−ut−1‖≤2C0α2log(1/α)+α2η.normsubscriptsubscript12subscript0superscript21superscript2\|u_t-u_t-1\|≤ 2C_0α^2 (1/α)+α^2η.∥ uitalic_t - uitalic_t - 1 ∥ ≤ 2 C0 α2 log ( 1 / α ) + α2 η . In total, we say that for t=t0,t0+1subscript0subscript01t=t_0,t_0+1t = t0 , t0 + 1, ‖ut+αHλ+κ¯(θt)−1∇L(θt)‖2+‖ut−ut−1‖21/2≤C3α2log(1α)superscriptsuperscriptnormsubscriptsubscript¯superscriptsubscript1∇subscript2superscriptnormsubscriptsubscript1212subscript3superscript21 \\|u_t+α H_λ+ κ( _t)^-1∇ L(% _t)\|^2+\|u_t-u_t-1\|^2 \^1/2≤ C_3α^2 % ( 1α ) ∥ uitalic_t + α Hitalic_λ + over¯ start_ARG κ end_ARG ( θitalic_t )- 1 ∇ L ( θitalic_t ) ∥2 + ∥ uitalic_t - uitalic_t - 1 ∥2 1 / 2 ≤ C3 α2 log ( divide start_ARG 1 end_ARG start_ARG α end_ARG ) for C3=3C1subscript33subscript1C_3=3C_1C3 = 3 C1 and sufficiently small α,ηα,ηα , η. Furthermore, by Lemma D.3 and Eq. (15) we have for all t∈[t0,2t0]subscript02subscript0t∈[t_0,2t_0]t ∈ [ t0 , 2 t0 ], ‖θt−θt−1‖≤C3ηα,normsubscriptsubscript1subscript3\| _t- _t-1\|≤ C_3ηα,∥ θitalic_t - θitalic_t - 1 ∥ ≤ C3 η α , as long as C3≥3K+5(1−μ)−1subscript335superscript11C_3≥ 3K+5(1-μ)^-1C3 ≥ 3 K + 5 ( 1 - μ )- 1. We are going to show by induction that maxt0≤t′≤t(∥ut′+αNλ¯(θt′)∥2+∥ut′−ut′−1∥2)1/2 _t_0≤ t ≤ t (\|u_t +α N_% λ( _t )\|^2+\|u_t -u_t % -1\|^2 )^1/2maxitalic_t start_POSTSUBSCRIPT 0 ≤ t′ ≤ t end_POSTSUBSCRIPT ( ∥ uitalic_t′ + α Nover¯ start_ARG λ end_ARG ( θitalic_t′ ) ∥2 + ∥ uitalic_t′ - uitalic_t′ - 1 ∥2 )1 / 2 ≤C3α2log(1α),absentsubscript3superscript21 ≤ C_3α^2 ( 1α ),≤ C3 α2 log ( divide start_ARG 1 end_ARG start_ARG α end_ARG ) , (17) maxt0≤t′≤t‖θt′−θt′−1‖subscriptsubscript0superscript′normsubscriptsuperscript′subscriptsuperscript′1 _t_0≤ t ≤ t\| _t - _t^% -1\|maxitalic_t start_POSTSUBSCRIPT 0 ≤ t′ ≤ t end_POSTSUBSCRIPT ∥ θitalic_t′ - θitalic_t′ - 1 ∥ ≤C3ηα,absentsubscript3 ≤ C_3ηα,≤ C3 η α , (18) We know that it holds at t=t0+1subscript01t=t_0+1t = t0 + 1. Given it holds for t, we need to extend it to t+11t+1t + 1. Then from Eq. 14, ut+1=(μ+1)ut−μut−1−η¯[α∇L(θt)+H(θt)ut]+O≤(9K3ηα2)subscript11subscriptsubscript1¯delimited-[]∇subscriptsubscriptsubscriptsubscript9superscript3superscript2u_t+1=(μ+1)u_t-μ u_t-1- η [α∇ L( _t% )+H( _t)u_t ]+O_≤(9K^3ηα^2)uitalic_t + 1 = ( μ + 1 ) uitalic_t - μ uitalic_t - 1 - over¯ start_ARG η end_ARG [ α ∇ L ( θitalic_t ) + H ( θitalic_t ) uitalic_t ] + O≤ ( 9 K3 η α2 ) Using the fact that ∇L(θ)=Hλ¯(θt)Nλ¯(θt)∇subscript¯subscriptsubscript¯subscript∇ L(θ)=H_ λ( _t)N_ λ(% _t)∇ L ( θ ) = Hover¯ start_ARG λ end_ARG ( θitalic_t ) Nover¯ start_ARG λ end_ARG ( θitalic_t ) by definition, we obtain the recursion ut+1+αNλ¯(θt)=(I−η¯H(θt)λ¯(θt))(ut+αNλ¯(θt))+μ(ut−ut−1)+O≤(9K3ηα2),subscript1subscript¯subscript¯subscriptsubscript¯subscriptsubscriptsubscript¯subscriptsubscriptsubscript1subscript9superscript3superscript2u_t+1+α N_ λ( _t)=(I- ηH( _% t)_ λ( _t))(u_t+α N_ λ(% _t))+μ(u_t-u_t-1)+O_≤(9K^3ηα^2),uitalic_t + 1 + α Nover¯ start_ARG λ end_ARG ( θitalic_t ) = ( I - over¯ start_ARG η end_ARG H ( θitalic_t )over¯ start_ARG λ end_ARG ( θitalic_t ) ) ( uitalic_t + α Nover¯ start_ARG λ end_ARG ( θitalic_t ) ) + μ ( uitalic_t - uitalic_t - 1 ) + O≤ ( 9 K3 η α2 ) , which can also be rewritten in the form (ut+1+αNλ¯(θt)ut+1−ut)=(I−η¯Hλ¯(θt)μ−η¯Hλ¯(θt)μ)(ut+αNλ¯(θt)ut−ut−1)+O≤(9K3ηα2)matrixsubscript1subscript¯subscriptsubscript1subscriptmatrix¯subscript¯subscript¯subscript¯subscriptmatrixsubscriptsubscript¯subscriptsubscriptsubscript1subscript9superscript3superscript2 pmatrixu_t+1+α N_ λ( _t)\\ u_t+1-u_t pmatrix= pmatrixI- ηH_ % λ( _t)&μ\\ - ηH_ λ( _t)&μ pmatrix % pmatrixu_t+α N_ λ( _t)\\ u_t-u_t-1 pmatrix+O_≤(9K^3ηα^2)( start_ARG start_ROW start_CELL uitalic_t + 1 + α Nover¯ start_ARG λ end_ARG ( θitalic_t ) end_CELL end_ROW start_ROW start_CELL uitalic_t + 1 - uitalic_t end_CELL end_ROW end_ARG ) = ( start_ARG start_ROW start_CELL I - over¯ start_ARG η end_ARG Hover¯ start_ARG λ end_ARG ( θitalic_t ) end_CELL start_CELL μ end_CELL end_ROW start_ROW start_CELL - over¯ start_ARG η end_ARG Hover¯ start_ARG λ end_ARG ( θitalic_t ) end_CELL start_CELL μ end_CELL end_ROW end_ARG ) ( start_ARG start_ROW start_CELL uitalic_t + α Nover¯ start_ARG λ end_ARG ( θitalic_t ) end_CELL end_ROW start_ROW start_CELL uitalic_t - uitalic_t - 1 end_CELL end_ROW end_ARG ) + O≤ ( 9 K3 η α2 ) From the proof of Lemma D.1, we know that the norm of the matrix in the middle is bounded by 1−1−μλη111- 1-μλη1 - divide start_ARG 1 - μ end_ARG start_ARG λ η end_ARG. Therefore, (‖ut+1+αNλ¯(θt)‖2+‖ut+1−ut‖2)1/2≤(1−λη1−μ)C3α2log(1α)+O≤(9K3ηα2)superscriptsuperscriptnormsubscript1subscript¯subscript2superscriptnormsubscript1subscript21211subscript3superscript21subscript9superscript3superscript2 (\|u_t+1+α N_ λ( _t)\|^2+\|u_t+1-u_t% \|^2 )^1/2≤ (1- λη1-μ )C_3α^2% ( 1α )+O_≤(9K^3ηα^2)( ∥ uitalic_t + 1 + α Nover¯ start_ARG λ end_ARG ( θitalic_t ) ∥2 + ∥ uitalic_t + 1 - uitalic_t ∥2 )1 / 2 ≤ ( 1 - divide start_ARG λ η end_ARG start_ARG 1 - μ end_ARG ) C3 α2 log ( divide start_ARG 1 end_ARG start_ARG α end_ARG ) + O≤ ( 9 K3 η α2 ) (19) To replace θtsubscript _tθitalic_t with θt+1subscript1 _t+1θitalic_t + 1 in the LHS, we need to bound the difference θt+1−θtsubscript1subscript _t+1- _tθitalic_t + 1 - θitalic_t, since by the regularity conditions and the triangle inequality we have (‖ut+1+αNλ¯(θt+1)‖2+‖ut+1−ut‖2)1/2≤(1−λη1−μ)C3α2log(1α)+9K3ηα2+Kα‖θt+1−θt‖superscriptsuperscriptdelimited-∥subscript1subscript¯subscript12superscriptdelimited-∥subscript1subscript21211subscript3superscript219superscript3superscript2delimited-∥subscript1subscript (\|u_t+1+α N_ λ( _t+1)\|^2+\|u_t+1-u_% t\|^2 )^1/2\\ ≤ (1- λη1-μ )C_3α^2 ( 1% α )+9K^3ηα^2+Kα\| _t+1- _t\|start_ROW start_CELL ( ∥ uitalic_t + 1 + α Nover¯ start_ARG λ end_ARG ( θitalic_t + 1 ) ∥2 + ∥ uitalic_t + 1 - uitalic_t ∥2 )1 / 2 end_CELL end_ROW start_ROW start_CELL ≤ ( 1 - divide start_ARG λ η end_ARG start_ARG 1 - μ end_ARG ) C3 α2 log ( divide start_ARG 1 end_ARG start_ARG α end_ARG ) + 9 K3 η α2 + K α ∥ θitalic_t + 1 - θitalic_t ∥ end_CELL end_ROW We first observe that from (19), the second term does not compensate −C3λη1−μα2log(1α)+9K3ηα2≤0subscript31superscript219superscript3superscript20-C_3 λη1-μα^2 ( 1α )+9K^% 3ηα^2≤ 0- C3 divide start_ARG λ η end_ARG start_ARG 1 - μ end_ARG α2 log ( divide start_ARG 1 end_ARG start_ARG α end_ARG ) + 9 K3 η α2 ≤ 0 (for small enough α), so that by the regularity condition, ‖ut+1‖≤Kα+C3α2log(1α)≤3Kαnormsubscript1subscript3superscript213\|u_t+1\|≤ Kα+C_3α^2 ( 1α )≤ 3Kα∥ uitalic_t + 1 ∥ ≤ K α + C3 α2 log ( divide start_ARG 1 end_ARG start_ARG α end_ARG ) ≤ 3 K α (20) We use (13) to bound θt+1−θtsubscript1subscript _t+1- _tθitalic_t + 1 - θitalic_t. Firstly, notice that given the approximation of utsubscriptu_tuitalic_t by Nλ¯(θt)subscript¯subscriptN_ λ( _t)Nover¯ start_ARG λ end_ARG ( θitalic_t ), we have that α∇L(θt)+Hλ(θt)ut=−(λ¯−λ)ut+α∇L(θt)+Hλ(θt)ut=O≤(3(1−μ)κ¯Kα+C3α2log(1/α))∇subscriptsubscriptsubscriptsubscript¯subscript∇subscriptsubscriptsubscriptsubscriptsubscript31¯subscript3superscript21α∇ L( _t)+H_λ( _t)u_t=-( λ-% λ)u_t+α∇ L( _t)+H_λ( _t)u_t=O_≤% (3(1-μ) κKα+C_3α^2 (1/α))α ∇ L ( θitalic_t ) + Hitalic_λ ( θitalic_t ) uitalic_t = - ( over¯ start_ARG λ end_ARG - λ ) uitalic_t + α ∇ L ( θitalic_t ) + Hitalic_λ ( θitalic_t ) uitalic_t = O≤ ( 3 ( 1 - μ ) over¯ start_ARG κ end_ARG K α + C3 α2 log ( 1 / α ) ) Then, we get from (13), ‖θt+1−θt‖≤μC3ηα+3(1−μ)κ¯Kα+C3α2log(1/α)+9K3ηα2≤C3ηα,normsubscript1subscriptsubscript331¯subscript3superscript219superscript3superscript2subscript3\| _t+1- _t\|≤μ C_3ηα+3(1-μ) κK% α+C_3α^2 (1/α)+9K^3ηα^2≤ C_3ηα,∥ θitalic_t + 1 - θitalic_t ∥ ≤ μ C3 η α + 3 ( 1 - μ ) over¯ start_ARG κ end_ARG K α + C3 α2 log ( 1 / α ) + 9 K3 η α2 ≤ C3 η α , where we assume that C3≥6κ¯Ksubscript36¯C_3≥ 6 κKC3 ≥ 6 over¯ start_ARG κ end_ARG K and α is sufficiently small. Hence, the induction assumption on the difference ‖θt+1−θt‖normsubscript1subscript\| _t+1- _t\|∥ θitalic_t + 1 - θitalic_t ∥ stands. We finalize the induction step by plugging this bound into the the inequality above, (‖ut+1+αNλ¯(θt+1)‖2+‖ut+1−ut‖2)1/2superscriptsuperscriptnormsubscript1subscript¯subscript12superscriptnormsubscript1subscript212 (\|u_t+1+α N_ λ( _t+1)\|^2% +\|u_t+1-u_t\|^2 )^1/2( ∥ uitalic_t + 1 + α Nover¯ start_ARG λ end_ARG ( θitalic_t + 1 ) ∥2 + ∥ uitalic_t + 1 - uitalic_t ∥2 )1 / 2 ≤(1−λη1−μ)C3α2log(1α)+9K3ηα2+C3Kηα2absent11subscript3superscript219superscript3superscript2subscript3superscript2 ≤ (1- λη1-μ )C_3α^2 % ( 1α )+9K^3ηα^2+C_3Kηα^2≤ ( 1 - divide start_ARG λ η end_ARG start_ARG 1 - μ end_ARG ) C3 α2 log ( divide start_ARG 1 end_ARG start_ARG α end_ARG ) + 9 K3 η α2 + C3 K η α2 =C3α2((1−λη1−μ)log(1/α)+(9K3+C3K)η)absentsubscript3superscript21119superscript3subscript3 =C_3α^2 ((1- λη1-μ) (1/α)% +(9K^3+C_3K)η )= C3 α2 ( ( 1 - divide start_ARG λ η end_ARG start_ARG 1 - μ end_ARG ) log ( 1 / α ) + ( 9 K3 + C3 K ) η ) ≤C3α2log(1α),absentsubscript3superscript21 ≤ C_3α^2 ( 1α ),≤ C3 α2 log ( divide start_ARG 1 end_ARG start_ARG α end_ARG ) , where the latter follows from the fact that λη1−μlog(1α)≥(9K3+C3K)η119superscript3subscript3 λη1-μ ( 1α )≥(9K^3+C_3K) start_ARG λ η end_ARG start_ARG 1 - μ end_ARG log ( divide start_ARG 1 end_ARG start_ARG α end_ARG ) ≥ ( 9 K3 + C3 K ) η for small enough α. Hence, we have shown by induction that (17), (18) hold for the rest of the training trajectory. Connecting to NG descent. Now we can finally show that updates mimic natural gradient descent, θ¯t+1=θ¯t−γHλ¯(θ¯t)−1∇L(θ¯t),θ¯0=θ0,formulae-sequencesubscript¯1subscript¯subscript¯superscriptsubscript¯1∇subscript¯subscript¯0subscript0 θ_t+1= θ_t-γ H_ λ(% θ_t)^-1∇ L( θ_t), % θ_0= _0,over¯ start_ARG θ end_ARGt + 1 = over¯ start_ARG θ end_ARGt - γ Hover¯ start_ARG λ end_ARG ( over¯ start_ARG θ end_ARGt )- 1 ∇ L ( over¯ start_ARG θ end_ARGt ) , over¯ start_ARG θ end_ARG0 = θ0 , where recall that γ=κ¯ηα¯γ= κηαγ = over¯ start_ARG κ end_ARG η α. Using Eq. (15) we have, θt+1−θ¯t+1subscript1subscript¯1 _t+1- θ_t+1θitalic_t + 1 - over¯ start_ARG θ end_ARGt + 1 =θt+κ¯ηut+(ut+1−ut)−θ¯t+κ¯ηαNλ¯(θ¯t)absentsubscript¯subscriptsubscript1subscriptsubscript¯subscript¯subscript¯ =θ_t+ κη u_t+(u_t+1-u_t)-% θ_t+ κηα N_ λ(% θ_t)= θitalic_t + over¯ start_ARG κ end_ARG η uitalic_t + ( uitalic_t + 1 - uitalic_t ) - over¯ start_ARG θ end_ARGt + over¯ start_ARG κ end_ARG η α Nover¯ start_ARG λ end_ARG ( over¯ start_ARG θ end_ARGt ) =θt−θ¯t+κ¯η(ut+αNλ¯(θ¯t))+(ut+1−ut)absentsubscriptsubscript¯subscriptsubscript¯subscript¯subscript1subscript =θ_t- θ_t+ κη(u_t+% α N_ λ( θ_t))+(u_t+1-u_t)= θitalic_t - over¯ start_ARG θ end_ARGt + over¯ start_ARG κ end_ARG η ( uitalic_t + α Nover¯ start_ARG λ end_ARG ( over¯ start_ARG θ end_ARGt ) ) + ( uitalic_t + 1 - uitalic_t ) =θt−θ¯t+κ¯η(ut+αNλ¯(θt))+κ¯ηα(Nλ¯(θ¯t)−Nλ¯(θt))+(ut+1−ut)absentsubscriptsubscript¯subscriptsubscript¯subscript¯subscript¯subscript¯subscript¯subscriptsubscript1subscript =θ_t- θ_t+ κη(u_t+% α N_ λ(θ_t))+ κηα(N_% λ( θ_t)-N_ λ(θ_t% ))+(u_t+1-u_t)= θitalic_t - over¯ start_ARG θ end_ARGt + over¯ start_ARG κ end_ARG η ( uitalic_t + α Nover¯ start_ARG λ end_ARG ( θitalic_t ) ) + over¯ start_ARG κ end_ARG η α ( Nover¯ start_ARG λ end_ARG ( over¯ start_ARG θ end_ARGt ) - Nover¯ start_ARG λ end_ARG ( θitalic_t ) ) + ( uitalic_t + 1 - uitalic_t ) Denoting Δt=‖θt−θ¯t−ut‖subscriptΔnormsubscriptsubscript¯subscript _t=\| _t- θ_t-u_t\|Δitalic_t = ∥ θitalic_t - over¯ start_ARG θ end_ARGt - uitalic_t ∥, we first obtain by the regularity condition and Eq. (20), ‖Nλ¯(θ¯t)−Nλ¯(θt)‖≤K‖θt−θ¯t‖≤K(Δt+3Kα).normsubscript¯subscript¯subscript¯subscriptnormsubscriptsubscript¯subscriptΔ3\|N_ λ( θ_t)-N_ λ(θ% _t)\|≤ K\|θ_t- θ_t\|≤ K( _t+3Kα).∥ Nover¯ start_ARG λ end_ARG ( over¯ start_ARG θ end_ARGt ) - Nover¯ start_ARG λ end_ARG ( θitalic_t ) ∥ ≤ K ∥ θitalic_t - over¯ start_ARG θ end_ARGt ∥ ≤ K ( Δitalic_t + 3 K α ) . Therefore, we get the recursive inequality, Δt+1subscriptΔ1 _t+1Δitalic_t + 1 ≤Δt+C3κ¯ηα2log1α+Kκ¯ηα(Δt+3Kα)absentsubscriptΔsubscript3¯superscript21¯subscriptΔ3 ≤ _t+C_3 κηα^2 1% α+K κηα( _t+3Kα)≤ Δitalic_t + C3 over¯ start_ARG κ end_ARG η α2 log divide start_ARG 1 end_ARG start_ARG α end_ARG + K over¯ start_ARG κ end_ARG η α ( Δitalic_t + 3 K α ) ≤(1+Kκ¯ηα)Δt+(C3log(1/α)+3K2)κ¯ηα2absent1¯subscriptΔsubscript313superscript2¯superscript2 ≤(1+K κηα) _t+(C_3 (1/% α)+3K^2) κηα^2≤ ( 1 + K over¯ start_ARG κ end_ARG η α ) Δitalic_t + ( C3 log ( 1 / α ) + 3 K2 ) over¯ start_ARG κ end_ARG η α2 ≤…absent… ≤…≤ … ≤(1+Kκ¯ηα)t+1−t0Δt0+(1+Kκ¯ηα)t+1−t0(C3αlog(1/α)+3K2α)κ¯ηα(t+1−t0)absentsuperscript1¯1subscript0subscriptΔsubscript0superscript1¯1subscript0subscript313superscript2¯1subscript0 ≤(1+K κηα)^t+1-t_0 _t_0+(1% +K κηα)^t+1-t_0(C_3α (1/α)+3K^2% α) κηα(t+1-t_0)≤ ( 1 + K over¯ start_ARG κ end_ARG η α )t + 1 - t0 Δitalic_t start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT + ( 1 + K over¯ start_ARG κ end_ARG η α )t + 1 - t0 ( C3 α log ( 1 / α ) + 3 K2 α ) over¯ start_ARG κ end_ARG η α ( t + 1 - t0 ) Given that the maximal number of steps T satisfies ηαT≤Aηα T≤ Aη α T ≤ A (which corresponds to the finite sum of the step sizes of NG, and it is one of the conditions of the theorem), we have that (1+Kκ¯ηα)t−t0≤(1+Kκ¯ηα)T≤exp(Kκ¯A),superscript1¯subscript0superscript1¯(1+K κηα)^t-t_0≤(1+K κηα)% ^T≤ (K κA),( 1 + K over¯ start_ARG κ end_ARG η α )t - t0 ≤ ( 1 + K over¯ start_ARG κ end_ARG η α )T ≤ exp ( K over¯ start_ARG κ end_ARG A ) , bounded by a constant. Therefore, for all t≤Tt≤ Tt ≤ T, Δt≤O(Δt0+αlog(1/α)).subscriptΔsubscriptΔsubscript01 _t≤ O( _t_0+α (1/α)).Δitalic_t ≤ O ( Δitalic_t start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT + α log ( 1 / α ) ) . It is left to bound Δt0subscriptΔsubscript0 _t_0Δitalic_t start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT, which simply relies on the fact that t0subscript0t_0t0 steps is not enough to diverge too far from the initial point θ0subscript0 _0θ0. Indeed, ‖θ¯t0−θ0‖≤γ∑t=0t0−1‖Nλ¯(θ¯t)‖≤Kt0γ=O(αlog(1/α)),normsubscript¯subscript0subscript0superscriptsubscript0subscript01normsubscript¯subscript¯subscript01\| θ_t_0- _0\|≤γ _t=0^t_0-1\|N_% λ( θ_t)\|≤ Kt_0γ=O(α (1/% α)),∥ over¯ start_ARG θ end_ARGt start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT - θ0 ∥ ≤ γ ∑t = 0t0 - 1 ∥ Nover¯ start_ARG λ end_ARG ( over¯ start_ARG θ end_ARGt ) ∥ ≤ K t0 γ = O ( α log ( 1 / α ) ) , and similarly, recalling Lemma D.3, ‖θ¯t0−θ0‖≤(t0−1)At≤t0(3Kκ¯ηα+5(1−μ¯)−1K2ηα)=O(αlog(1/α)).normsubscript¯subscript0subscript0subscript01subscriptsubscript03¯5superscript1¯1superscript21\| θ_t_0- _0\|≤(t_0-1)A_t≤ t_0(3K% κηα+5(1- μ)^-1K^2ηα)=O(α% (1/α)).∥ over¯ start_ARG θ end_ARGt start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT - θ0 ∥ ≤ ( t0 - 1 ) Aitalic_t ≤ t0 ( 3 K over¯ start_ARG κ end_ARG η α + 5 ( 1 - over¯ start_ARG μ end_ARG )- 1 K2 η α ) = O ( α log ( 1 / α ) ) . We also note that Δt0≤‖θ¯t0−θ0‖+‖θt0−θ0‖subscriptΔsubscript0normsubscript¯subscript0subscript0normsubscriptsubscript0subscript0 _t_0≤\| θ_t_0- _0\|+\|θ_t_0-% _0\|Δitalic_t start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT ≤ ∥ over¯ start_ARG θ end_ARGt start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT - θ0 ∥ + ∥ θitalic_t start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT - θ0 ∥. Then, ‖θt−θ¯t‖≤‖ut‖+Δt=O(α)+O(αlog(1/α))=O(αlog(1/α)).normsubscriptsubscript¯normsubscriptsubscriptΔ11\| _t- θ_t\|≤\|u_t\|+ _t=O(α)+O(% α (1/α))=O(α (1/α)).∥ θitalic_t - over¯ start_ARG θ end_ARGt ∥ ≤ ∥ uitalic_t ∥ + Δitalic_t = O ( α ) + O ( α log ( 1 / α ) ) = O ( α log ( 1 / α ) ) .