One recurring theme in empirical deep learning is that we want to keep intermediate activations at a healthy size – we don’t want them to blow up or vanish.

This principle has inspired a bunch of design choices that are popular in empirical deep learning today – for example initialization schemes (e.g., Xavier and Kaiming He initialization), layer normalization, and gradient clipping.

Muon[1]Muon: An optimizer for hidden layers in neural networks [link]
Keller Jordan, Yuchen Jin, Vlado Boza, Jiacheng You, Franz Cesista, Laker Newhouse, and Jeremy Bernstein. 2024.
is inspired by this same idea, applied to optimization. The goal is to keep activations at every layer reasonably sized (we’ll define “reasonably sized” in a moment); these activations should be reasonably sized not just at initialization, but also throughout training.

Natural norms

We first need to define what we mean when we say that an activation vector v∈Rd\mathbf{v} \in \mathbb{R}^d is “reasonably sized”.

One standard way to measure the size of a vector is the ℓ2\ell_2 norm: ∥v∥2=∑i=1dvi2\|\mathbf{v}\|_2 = \sqrt{\sum_{i=1}^d v_i^2}.

However, this norm is dimension-dependent: for example, a vector with all entries ±1\pm 1 has ∥v∥2=d\|\mathbf{v}\|_2 = \sqrt{d}, which grows with dimension.

The natural norm for vectors: RMS norm

The root mean square (RMS) norm fixes this by averaging over entries instead of summing:

∥v∥RMS:=1d∑i=1dvi2=1d∥v∥2.\|\mathbf{v}\|_{\text{RMS}} := \sqrt{\frac{1}{d} \sum_{i=1}^d v_i^2} = \frac{1}{\sqrt{d}} \|\mathbf{v}\|_2.

Now a ±1\pm 1 vector has ∥v∥RMS=1\|\mathbf{v}\|_{\text{RMS}} = 1 regardless of dimension.

Sanity check with a Gaussian random vector
E[∥v∥RMS2]=E ⁣[(1d∥v∥2)2]=1dE[∥v∥22]=1dE ⁣[∑i=1dvi2]=1d∑i=1dE[vi2]=1d⋅d⋅1=1,\begin{aligned} \mathbb{E}[\|\mathbf{v}\|_{\text{RMS}}^2] &= \mathbb{E}\!\left[\left(\frac{1}{\sqrt{d}} \|\mathbf{v}\|_2\right)^2\right] \\ &= \frac{1}{d}\mathbb{E}[\|\mathbf{v}\|_2^2] \\ &= \frac{1}{d}\mathbb{E}\!\left[\sum_{i=1}^d v_i^2\right] \\ &= \frac{1}{d}\sum_{i=1}^d \mathbb{E}[v_i^2] \\ &= \frac{1}{d} \cdot d \cdot 1 = 1, \end{aligned}

since each vi∼N(0,1)v_i \sim \mathcal{N}(0, 1) has E[vi2]=Var(vi)=1\mathbb{E}[v_i^2] = \text{Var}(v_i) = 1.

Compare with the Euclidean norm:

E[∥v∥22]=E ⁣[∑i=1dvi2]=∑i=1dE[vi2]=d⋅1=d,\mathbb{E}[\|\mathbf{v}\|_2^2] = \mathbb{E}\!\left[\sum_{i=1}^d v_i^2\right] = \sum_{i=1}^d \mathbb{E}[v_i^2] = d \cdot 1 = d,

which grows with dimension.

Intuitively, the RMS norm measures the scale of a typical entry of v\mathbf{v}, rather than accumulating across all entries like the ℓ2\ell_2 norm does. This makes it the natural norm for dense vectors1: ∥v∥RMS=Θ(1)\|\mathbf{v}\|_{\text{RMS}} = \Theta(1) exactly captures “the entries of v\mathbf{v} are Θ(1)\Theta(1)-sized”, independent of dimension.

The natural norm for matrices: RMS-to-RMS operator norm

What does it mean for a weight matrix W∈Rdout×din\mathbf{W} \in \mathbb{R}^{d_\text{out} \times d_\text{in}} to be “well-behaved”? It should map reasonably-sized inputs to reasonably-sized outputs: if ∥x∥RMS=Θ(1)\|\mathbf{x}\|_{\text{RMS}} = \Theta(1), then we want ∥Wx∥RMS=Θ(1)\|\mathbf{W}\mathbf{x}\|_{\text{RMS}} = \Theta(1).

The natural way to measure this is the RMS-to-RMS operator norm – the worst-case RMS stretch:

∥W∥RMS→RMS:=max⁡x≠0∥Wx∥RMS∥x∥RMS.\|\mathbf{W}\|_{\text{RMS} \to \text{RMS}} := \max_{\mathbf{x} \neq \mathbf{0}} \frac{\|\mathbf{W}\mathbf{x}\|_{\text{RMS}}}{\|\mathbf{x}\|_{\text{RMS}}}.

This is the natural spectral norm. It relates to the standard spectral norm ∥W∥∗\|\mathbf{W}\|_* (largest singular value) by a dimensional factor:

∥W∥RMS→RMS=dindout⋅∥W∥∗.\|\mathbf{W}\|_{\text{RMS} \to \text{RMS}} = \sqrt{\frac{d_\text{in}}{d_\text{out}}} \cdot \|\mathbf{W}\|_*.
Relating RMS-to-RMS and spectral norm
∥Wx∥RMS∥x∥RMS=1dout∥Wx∥21din∥x∥2=dindout⋅∥Wx∥2∥x∥2.\frac{\|\mathbf{W}\mathbf{x}\|_{\text{RMS}}}{\|\mathbf{x}\|_{\text{RMS}}} = \frac{\frac{1}{\sqrt{d_\text{out}}} \|\mathbf{W}\mathbf{x}\|_2}{\frac{1}{\sqrt{d_\text{in}}} \|\mathbf{x}\|_2} = \sqrt{\frac{d_\text{in}}{d_\text{out}}} \cdot \frac{\|\mathbf{W}\mathbf{x}\|_2}{\|\mathbf{x}\|_2}.

Taking the max over x\mathbf{x}:

∥W∥RMS→RMS=dindout⋅∥W∥∗.\|\mathbf{W}\|_{\text{RMS} \to \text{RMS}} = \sqrt{\frac{d_\text{in}}{d_\text{out}}} \cdot \|\mathbf{W}\|_*.

The spectral scaling condition

Let’s recall what we’re trying to achieve. We want every layer’s weight matrix W\mathbf{W} to map reasonably-sized inputs to reasonably-sized outputs – that is, ∥W∥RMS→RMS=Θ(1)\|\mathbf{W}\|_{\text{RMS} \to \text{RMS}} = \Theta(1). Using the relationship we derived above, this is equivalent to:

∥W∥∗=Θ ⁣(doutdin).\|\mathbf{W}\|_* = \Theta\!\left(\sqrt{\frac{d_\text{out}}{d_\text{in}}}\right).

We want this to hold not just at initialization, but throughout training. After a weight update W→W+ΔW\mathbf{W} \to \mathbf{W} + \Delta \mathbf{W}, the change in output for input x\mathbf{x} is Δy=ΔW x\Delta \mathbf{y} = \Delta \mathbf{W} \, \mathbf{x}. By definition of the operator norm:

∥Δy∥RMS=∥ΔW x∥RMS≤∥ΔW∥RMS→RMS⋅∥x∥RMS.\|\Delta \mathbf{y}\|_{\text{RMS}} = \|\Delta \mathbf{W} \, \mathbf{x}\|_{\text{RMS}} \leq \|\Delta \mathbf{W}\|_{\text{RMS} \to \text{RMS}} \cdot \|\mathbf{x}\|_{\text{RMS}}.

So if the input is reasonably sized (∥x∥RMS=Θ(1)\|\mathbf{x}\|_{\text{RMS}} = \Theta(1)) and we want the output change to also be reasonably sized (∥Δy∥RMS=Θ(1)\|\Delta \mathbf{y}\|_{\text{RMS}} = \Theta(1)), we need ∥ΔW∥RMS→RMS=Θ(1)\|\Delta \mathbf{W}\|_{\text{RMS} \to \text{RMS}} = \Theta(1), i.e.:

∥ΔW∥∗=Θ ⁣(doutdin).\|\Delta \mathbf{W}\|_* = \Theta\!\left(\sqrt{\frac{d_\text{out}}{d_\text{in}}}\right).

Together, these two requirements form the spectral scaling condition from Yang et al. [2]A spectral condition for feature learning [link]
Greg Yang, James B. Simon, and Jeremy Bernstein. arXiv preprint. 2024.
: both the weights and their updates should have Θ(1)\Theta(1) RMS-to-RMS operator norm, or equivalently, spectral norm Θ(dout/din)\Theta(\sqrt{d_\text{out}/d_\text{in}}).

Deriving Muon

We’ve established that ΔW\Delta \mathbf{W} should have Θ(1)\Theta(1) RMS-to-RMS operator norm. That is, ∥ΔW∥RMS→RMS≤α\|\Delta \mathbf{W}\|_{\text{RMS} \to \text{RMS}} \leq \alpha for some constant α\alpha. Since ∥W∥RMS→RMS=din/dout⋅∥W∥∗\|\mathbf{W}\|_{\text{RMS} \to \text{RMS}} = \sqrt{d_\text{in}/d_\text{out}} \cdot \|\mathbf{W}\|_*, this is equivalent to ∥ΔW∥∗≤η\|\Delta \mathbf{W}\|_* \leq \eta where η:=αdout/din\eta := \alpha \sqrt{d_\text{out}/d_\text{in}}.

Given this budget, we want to decrease the loss as much as possible. To first order, the change in loss is ⟨∇WL,ΔW⟩\langle \nabla_{\mathbf{W}} \mathcal{L}, \Delta \mathbf{W} \rangle, where ⟨A,B⟩:=∑ijAijBij\langle \mathbf{A}, \mathbf{B} \rangle := \sum_{ij} A_{ij} B_{ij} is the entrywise inner product between matrices. So we want to solve:

min⁡ΔW  ⟨∇WL,ΔW⟩subject to∥ΔW∥∗≤η.\min_{\Delta \mathbf{W}} \; \langle \nabla_{\mathbf{W}} \mathcal{L}, \Delta \mathbf{W} \rangle \quad \text{subject to} \quad \|\Delta \mathbf{W}\|_* \leq \eta.

We can write the gradient in its SVD:

∇WL=∑i=1rσiuivi⊤,\nabla_{\mathbf{W}} \mathcal{L} = \sum_{i=1}^r \sigma_i \mathbf{u}_i \mathbf{v}_i^\top,

with singular values σ1≥σ2≥⋯≥σr>0\sigma_1 \geq \sigma_2 \geq \cdots \geq \sigma_r > 0.

Plugging in the SVD, and using linearity of inner product, we get:

⟨∇WL,ΔW⟩  =  ⟨∑i=1rσiuivi⊤,  ΔW⟩  =  ∑i=1rσi ⟨uivi⊤,ΔW⟩.\langle \nabla_{\mathbf{W}} \mathcal{L}, \Delta \mathbf{W} \rangle \;=\; \left\langle \sum_{i=1}^r \sigma_i \mathbf{u}_i \mathbf{v}_i^\top, \; \Delta \mathbf{W} \right\rangle \;=\; \sum_{i=1}^r \sigma_i \, \langle \mathbf{u}_i \mathbf{v}_i^\top, \Delta \mathbf{W} \rangle.

Note that only the projection of ΔW\Delta \mathbf{W} onto the singular directions {uivi⊤}\{\mathbf{u}_i \mathbf{v}_i^\top\} appears in the objective – any orthogonal component doesn’t help decrease the loss, but can only increase ∥ΔW∥∗\|\Delta \mathbf{W}\|_*. Thus, without loss of generality, we can write ΔW=∑i=1rciuivi⊤\Delta \mathbf{W} = \sum_{i=1}^r c_i \mathbf{u}_i \mathbf{v}_i^\top.

This is an SVD of ΔW\Delta \mathbf{W} with singular values ∣ci∣\lvert c_i \rvert, so ∥ΔW∥∗=max⁡i∣ci∣\|\Delta \mathbf{W}\|_* = \max_i \lvert c_i \rvert, and the optimization reduces to:

min⁡c1,…,cr  ∑i=1rσi cisubject tomax⁡i∣ci∣≤η.\min_{c_1, \ldots, c_r} \; \sum_{i=1}^r \sigma_i \, c_i \quad \text{subject to} \quad \max_i \lvert c_i \rvert \leq \eta.

Since all σi>0\sigma_i > 0, the minimum is achieved by setting ci=−ηc_i = -\eta for all ii, giving:

ΔW=−η∑i=1ruivi⊤=−η UV⊤.\Delta \mathbf{W} = -\eta \sum_{i=1}^r \mathbf{u}_i \mathbf{v}_i^\top = -\eta \, \mathbf{U}\mathbf{V}^\top.

This is the Muon update: take the gradient’s SVD, replace all singular values with 1, and then scale by −η-\eta.

Comparison with gradient descent

What if we had constrained the Frobenius norm instead of the spectral norm? That is, what if we solved:

min⁡ΔW  ⟨∇WL,ΔW⟩subject to∥ΔW∥F≤η.\min_{\Delta \mathbf{W}} \; \langle \nabla_{\mathbf{W}} \mathcal{L}, \Delta \mathbf{W} \rangle \quad \text{subject to} \quad \|\Delta \mathbf{W}\|_F \leq \eta.

If we treat ΔW\Delta \mathbf{W} as a flat vector of entries, the Frobenius norm is just the ℓ2\ell_2 norm, and this is asking: which direction decreases the loss most per unit ℓ2\ell_2 step? That’s the definition of the gradient. So the solution is ΔW∝−∇WL\Delta \mathbf{W} \propto -\nabla_{\mathbf{W}} \mathcal{L}: gradient descent.

So the key difference between gradient descent and Muon stems from the norm constraint. Gradient descent effectively constrains the Frobenius norm, which treats the weight matrix as an unstructured vector of numbers; Muon effectively constrains the spectral norm, which measures how the matrix acts on inputs.

The spectral norm only constrains the largest singular value of ΔW\Delta \mathbf{W} – once you’ve spent your budget η\eta on the top singular direction, the remaining directions are free. Muon takes advantage of this by stepping equally in every singular direction. Gradient descent can’t do this: under the Frobenius norm, every direction draws from a shared budget (∑ci2≤η2\sum c_i^2 \leq \eta^2), so stepping more in one direction means stepping less in another. The optimal allocation sets ci∝−σic_i \propto -\sigma_i, concentrating on the largest singular directions at the expense of the smaller ones.

Making it practical: Newton-Schulz

The (ideal) Muon update requires computing UV⊤\mathbf{U}\mathbf{V}^\top from the gradient ∇WL=UΣV⊤\nabla_{\mathbf{W}} \mathcal{L} = \mathbf{U}\boldsymbol{\Sigma} \mathbf{V}^\top. Computing the full SVD at every step is too expensive. But we don’t need the full SVD – we just need the map UΣV⊤↦UV⊤\mathbf{U}\boldsymbol{\Sigma} \mathbf{V}^\top \mapsto \mathbf{U}\mathbf{V}^\top.

Turns out there are some algorithms that do exactly this – Newton-Schulz iterations. The algorithm works iteratively: at each step, it applies an operation to the matrix that acts independently on each singular value while preserving the singular vectors, and the operation is chosen so that all singular values converge to 1. In practice, 5-10 iterations suffice, and each iteration is just a few matrix multiplications – fast on GPUs. See these blog posts[3, 1]Deriving Muon [link]
Jeremy Bernstein. 2025.
Muon: An optimizer for hidden layers in neural networks [link]
Keller Jordan, Yuchen Jin, Vlado Boza, Jiacheng You, Franz Cesista, Laker Newhouse, and Jeremy Bernstein. 2024.
for more details.

The Muon update rule

Putting it all together, the Muon update for a weight matrix W∈Rdout×din\mathbf{W} \in \mathbb{R}^{d_\text{out} \times d_\text{in}} is:

W←W−α⋅doutdin⋅NewtonSchulz(∇WL).\mathbf{W} \leftarrow \mathbf{W} - \alpha \cdot \sqrt{\frac{d_\text{out}}{d_\text{in}}} \cdot \text{NewtonSchulz}(\nabla_{\mathbf{W}} \mathcal{L}).

Let’s unpack this:

  • NewtonSchulz(∇WL)\text{NewtonSchulz}(\nabla_{\mathbf{W}} \mathcal{L}) approximately orthogonalizes the gradient: UΣV⊤↦UV⊤\mathbf{U}\boldsymbol{\Sigma} \mathbf{V}^\top \mapsto \mathbf{U}\mathbf{V}^\top. This has spectral norm 11.
  • The dout/din\sqrt{d_\text{out}/d_\text{in}} factor scales the RMS-to-RMS learning rate α\alpha to the spectral norm budget: recall η=αdout/din\eta = \alpha \sqrt{d_\text{out}/d_\text{in}}.

Note that the dout/din\sqrt{d_\text{out}/d_\text{in}} factor absorbs the layer dimensions into the update, which helps α\alpha transfer across weight matrices of different shapes and sizes.

References

References are listed in order of citation.

  1. Muon: An optimizer for hidden layers in neural networks [link]
    Keller Jordan, Yuchen Jin, Vlado Boza, Jiacheng You, Franz Cesista, Laker Newhouse, and Jeremy Bernstein. 2024.
  2. A spectral condition for feature learning [link]
    Greg Yang, James B. Simon, and Jeremy Bernstein. arXiv preprint. 2024.
  3. Deriving Muon [link]
    Jeremy Bernstein. 2025.

Footnotes

  1. A dense vector is one where every entry contributes a comparable amount to the squared norm. By contrast, a sparse vector has only a constant number of non-negligible entries, regardless of dimension – e.g., a one-hot encoding vector. For sparse vectors, the ordinary ℓ2\ell_2 norm is already dimension-independent (e.g., a one-hot vector has ∥ei∥2=1\|\mathbf{e}_i\|_2 = 1, regardless of dimension), so no correction is needed. The RMS norm is the natural norm specifically for dense vectors. Throughout this note we focus on dense vectors, since hidden activations in a transformer can be thought of as dense; the embedding layer (which maps sparse inputs to dense activations) requires separate treatment. ↩