and
$\mathbf{w}^{*}_{S}=\arg\min_{\mathbf{w}\in\mathcal{W}}F_{S}(\mathbf{w},\bar{\mathbf{v}}_{T})$.
Similar to (C.2), we can show

$\displaystyle\frac{1}{T}\sum_{j=1}^{T}F_{S}(\mathbf{w}_{j},\mathbf{v}_{j})-F_{S}(\mathbf{w}^{*}_{S},\bar{\mathbf{v}}_{T})$
$\displaystyle\leq\frac{3G^{2}\eta}{T}\sum_{j=1}^{T}k_{j}+\frac{1}{T}\sum_{j=1}^{T}\big{(}F_{S}(\mathbf{w}_{j-k_{j}},\mathbf{v}_{j-k_{j}})-F_{S}(\mathbf{w}^{*}_{S},\mathbf{v}_{j-k_{j}})\big{)}.$
(C.17)

Let
$\xi_{j}=f(\mathbf{w}_{j-k_{j}},\mathbf{v}_{j-k_{j}};z_{i_{j}})-f(\mathbf{w}_{S}^{*},\mathbf{v}_{j-k_{j}};z_{i_{j}}).$
Observe that $|\xi_{j}-\mathbb{E}_{i_{j}}[\xi_{j}]|\leq 2B$. Then, applying
Lemma A.5 implies, with probability at least $1-\gamma/2$, that

$\sum_{j=1}^{T}\mathbb{E}_{i_{j}}[\xi_{j}]-\sum_{j=1}^{T}\xi_{j}\leq
2B\sqrt{2T\log(2/\gamma)}.$ (C.18)

Combining (C.2) and (C.18) together, we get

$\displaystyle\sum_{j=1}^{T}[F_{S}(\mathbf{w}_{j-k_{j}},\mathbf{v}_{j-k_{j}})\\!-\\!F_{S}(\mathbf{w}^{*}_{S},\mathbf{v}_{j-k_{j}})]\\!+\\!\sum_{j=1}^{T}\sum_{i=1}^{n}\big{(}[P^{k_{j}}]_{i_{j-k_{j}},i}\\!-\\!\frac{1}{n}\big{)}[f(\mathbf{w}_{j-k_{j}},\mathbf{v}_{j-k_{j}};z_{i})\\!-\\!f(\mathbf{w}_{S}^{*},\mathbf{v}_{j-k_{j}};z_{i})]$
$\displaystyle=\sum_{j=1}^{T}\mathbb{E}_{i_{j}}[f(\mathbf{w}_{j-k_{j}},\mathbf{v}_{j-k_{j}};z_{i_{j}})\\!-\\!f(\mathbf{w}^{*}_{S},\mathbf{v}_{j-k_{j}};z_{i_{j}})|\\{\mathbf{w}_{0},\mathbf{v}_{0}\\},\ldots,\\{\mathbf{w}_{j-k_{j}},\mathbf{v}_{j-k_{j}}\\},z_{i_{1}},\ldots,z_{i_{j-k_{j}}}]$
$\displaystyle\leq\sum_{j=1}^{T}[f(\mathbf{w}_{j-k_{j}},\mathbf{v}_{j-k_{j}};z_{i_{j}})-f(\mathbf{w}^{*}_{S},\mathbf{v}_{j-k_{j}};z_{i_{j}})]+2B\sqrt{2T\log(2/\gamma)}$

with probability at least $1-\gamma/2$. Putting (C.11) and (C.2) back into the
above inequality, we obtain

$\displaystyle\sum_{j=1}^{T}[F_{S}(\mathbf{w}_{j-k_{j}},\mathbf{v}_{j-k_{j}})-F_{S}(\mathbf{w}^{*}_{S},\mathbf{v}_{j-k_{j}})]$
$\displaystyle\leq\sum_{j=1}^{T}\sum_{i=1}^{n}\Big{(}\frac{1}{n}-[P^{k_{j}}]_{i_{j-k_{j}},i}\Big{)}[f(\mathbf{w}_{j-k_{j}},\mathbf{v}_{j-k_{j}};z_{i})-f(\mathbf{w}^{*}_{S},\mathbf{v}_{j-k_{j}};z_{i})]$
$\displaystyle\quad+\sum_{j=1}^{t}[f(\mathbf{w}_{j-k_{j}},\mathbf{v}_{j-k_{j}};z_{i_{j}})-f(\mathbf{w}^{*}_{S},\mathbf{v}_{j-k_{j}};z_{i_{j}})]+2B\sqrt{2T\log(2/\gamma)}$
$\displaystyle\leq\frac{D_{\mathbf{w}}^{2}+6G^{2}\eta^{2}\sum_{j=1}^{T}k_{j}+TG^{2}\eta^{2}}{2\eta}+2GK_{P}D_{\mathbf{w}}+G\sum_{j=K_{P}}^{T}\frac{1}{2j^{2}}+2B\sqrt{2T\log(2/\gamma)}.$
(C.19)

Now, plugging (C.2) back into (C.17), with probability at least $1-\gamma/2$,
there holds

$\displaystyle\frac{1}{T}\sum_{j=1}^{T}F_{S}(\mathbf{w}_{j},\mathbf{v}_{j})-\min_{\mathbf{w}\in\mathcal{W}}F_{S}(\mathbf{w},\bar{\mathbf{v}}_{T})$
$\displaystyle\leq\frac{6G^{2}\eta\sum_{j=1}^{T}k_{j}+2GK_{P}D_{\mathbf{w}}+G\sum_{j=K_{P}}^{T}\frac{1}{2j^{2}}}{T}+\frac{D_{\mathbf{w}}^{2}}{2T\eta}+\frac{G^{2}\eta}{2}+\frac{2B\sqrt{2\log(\frac{2}{\gamma})}}{\sqrt{T}}.$

In a similar way, we can show, with probability at least $1-\gamma/2$, that

$\displaystyle\max_{\mathbf{v}\in\mathcal{V}}F_{S}(\bar{\mathbf{w}}_{T},\mathbf{v})-\frac{1}{T}\sum_{j=1}^{T}F_{S}(\mathbf{w}_{j},\mathbf{v}_{j})$
$\displaystyle\leq\frac{6G^{2}\eta\sum_{j=1}^{T}k_{j}+2GK_{P}D_{\mathbf{v}}+G\sum_{j=K_{P}}^{T}\frac{1}{2j^{2}}}{T}+\frac{D_{\mathbf{v}}^{2}}{2T\eta}+\frac{G^{2}\eta}{2}+\frac{2B\sqrt{2\log(\frac{2}{\gamma})}}{\sqrt{T}}.$

Combining the above two inequalities together, with probability at least
$1-\gamma$, we get

$\displaystyle\max_{\mathbf{v}\in\mathcal{V}}F_{S}(\bar{\mathbf{w}}_{T},\mathbf{v})-\min_{\mathbf{w}\in\mathcal{W}}F_{S}(\mathbf{w},\bar{\mathbf{v}}_{T})$
$\displaystyle\leq\frac{12G^{2}\eta\sum_{j=1}^{T}k_{j}\\!+\\!2GK_{P}(D_{\mathbf{w}}\\!+\\!D_{\mathbf{v}})\\!+\\!G\sum_{j=K_{P}}^{T}\frac{1}{j^{2}}}{T}\\!+\\!\frac{(D_{\mathbf{w}}+D_{\mathbf{v}})^{2}}{2T\eta}\\!+\\!G^{2}\eta\\!+\\!\frac{4B\sqrt{2\log(\frac{2}{\gamma})}}{\sqrt{T}}.$

Further, if we select $\eta\asymp 1/{\sqrt{T\log(T)}}$, according to Eq.(C.15)
we have

$\displaystyle\max_{\mathbf{v}\in\mathcal{V}}F_{S}(\bar{\mathbf{w}}_{T},\mathbf{v})-\min_{\mathbf{w}\in\mathcal{W}}F_{S}(\mathbf{w},\bar{\mathbf{v}}_{T})$
$\displaystyle=\mathcal{O}\Big{(}\frac{K_{P}}{T}+\eta+\frac{1+\sum_{j=1}^{T}k_{j}\eta^{2}}{T\eta}+\frac{B\sqrt{\log(1/{\gamma})}}{\sqrt{T}}\Big{)}$
$\displaystyle=\mathcal{O}\Big{(}\frac{\sqrt{\log(T)}}{\sqrt{T}}\Big{(}\frac{1}{\log(1/\lambda(P))}+B\sqrt{\log(1/\gamma)}\Big{)}+\frac{K_{P}}{T\min\\{1,\sqrt{C_{P}n\lambda(P)^{K_{P}}T\log(T)}\\}}\Big{)}$
$\displaystyle=\mathcal{O}\Big{(}\frac{\sqrt{\log(T)}}{\sqrt{T}}\Big{(}\frac{1}{\log(1/\lambda(P))}+B\sqrt{\log(1/\gamma)}\Big{)}\Big{)},$

where in the last equality we used $K_{P}=0$ due to $P=P^{\top}$. This
completes the proof. ∎

###  C.3 Proofs of Theorem 9 and Theorem 10

###### Proof of Theorem 9.

We can choose $\eta$ such that $T\eta^{2}\leq 1/(2L^{2})$ and therefore
Theorem 7 applies. According to part (a) of Theorem 7 we have

$\triangle^{w}(\bar{\mathbf{w}}_{T},\bar{\mathbf{v}}_{T})-\triangle^{w}_{\text{emp}}(\bar{\mathbf{w}}_{T},\bar{\mathbf{v}}_{T})\leq\frac{4G^{2}\sqrt{T}\eta}{\sqrt{n}}+\frac{8\sqrt{2}G^{2}T\eta}{n}.$

Combining the above inequality with Theorem C.3 together, we get

$\displaystyle\triangle^{w}(\bar{\mathbf{w}}_{T},\bar{\mathbf{v}}_{T})=\triangle^{w}(\bar{\mathbf{w}}_{T},\bar{\mathbf{v}}_{T})-\triangle^{w}_{\text{emp}}(\bar{\mathbf{w}}_{T},\bar{\mathbf{v}}_{T})$
$\displaystyle\leq\frac{4G^{2}\sqrt{T}\eta}{\sqrt{n}}+\frac{8\sqrt{2}G^{2}T\eta}{n}+G^{2}\eta+\frac{D^{2}}{2T\eta}+\frac{2GK_{P}D+12G^{2}\eta\sum_{j=1}^{T}k_{j}+G\sum_{j=K_{P}}^{T}1/{j^{2}}}{T}$

where $D=D_{\mathbf{w}}+D_{\mathbf{v}}$. If we choose $T\asymp n$ and
$\eta\asymp(T\log(T))^{-\frac{1}{2}}$, according to (C.15) we get

$\triangle^{w}(\bar{\mathbf{w}}_{T},\bar{\mathbf{v}}_{T})=\mathcal{O}\Big{(}\frac{\log(n)}{\sqrt{n}\log(1/\lambda(P))}+\frac{K_{P}}{n\min\\{1,n\sqrt{\log(n)C_{P}\lambda(P)^{K_{P}}}\\}}\Big{)}.$

Note Assumption 4 implies $K_{P}=0$, the proof of part (a) is completed.

Part (b) can be proved in a similar way (e.g., by combining part (b) of
Theorem 7 and Theorem C.3 together). We omit the proof for brevity. ∎

###### Proof of Theorem 10.

We use the following decomposition

$\displaystyle R(\bar{\mathbf{w}}_{T})-R(\mathbf{w}^{*})=$
$\displaystyle\Big{(}R(\bar{\mathbf{w}}_{T})-R_{S}(\bar{\mathbf{w}}_{T})\Big{)}+\Big{(}R_{S}(\bar{\mathbf{w}}_{T})-F_{S}(\mathbf{w}^{*},\bar{\mathbf{v}}_{T})\Big{)}$
$\displaystyle+\Big{(}F_{S}(\mathbf{w}^{*},\bar{\mathbf{v}}_{T})-F(\mathbf{w}^{*},\bar{\mathbf{v}}_{T})\Big{)}+\Big{(}F(\mathbf{w}^{*},\bar{\mathbf{v}}_{T})-R(\mathbf{w}^{*})\Big{)}.$

Note that $F(\mathbf{w}^{*},\bar{\mathbf{v}}_{T})\leq
F(\mathbf{w}^{*},{\mathbf{v}}^{*})$. Then we have

$\displaystyle R(\bar{\mathbf{w}}_{T})-R(\mathbf{w}^{*})\leq$
$\displaystyle\Big{(}R(\bar{\mathbf{w}}_{T})-R_{S}(\bar{\mathbf{w}}_{T})\Big{)}+\Big{(}R_{S}(\bar{\mathbf{w}}_{T})-F_{S}(\mathbf{w}^{*},\bar{\mathbf{v}}_{T})\Big{)}$
$\displaystyle+\Big{(}F_{S}(\mathbf{w}^{*},\bar{\mathbf{v}}_{T})-F(\mathbf{w}^{*},\bar{\mathbf{v}}_{T})\Big{)}.$

Taking the expectation on both sides gives

$\displaystyle\mathbb{E}_{S,\mathcal{A}}[R(\bar{\mathbf{w}}_{T})-R(\mathbf{w}^{*})]\leq$
$\displaystyle\mathbb{E}_{S,\mathcal{A}}[R(\bar{\mathbf{w}}_{T})-R_{S}(\bar{\mathbf{w}}_{T})]+\mathbb{E}_{S,\mathcal{A}}[R_{S}(\bar{\mathbf{w}}_{T})-F_{S}(\mathbf{w}^{*},\bar{\mathbf{v}}_{T})]$
$\displaystyle+\mathbb{E}_{S,\mathcal{A}}[F_{S}(\mathbf{w}^{*},\bar{\mathbf{v}}_{T})-F(\mathbf{w}^{*},\bar{\mathbf{v}}_{T})].$
(C.20)

According to part (a) of Theorem 8 we know

$\mathbb{E}_{S,\mathcal{A}}[R(\bar{\mathbf{w}}_{T})-R_{S}(\bar{\mathbf{w}}_{T})]\leq
4G^{2}(1+L/\rho)\Big{(}\frac{\sqrt{T}\eta}{\sqrt{n}}+\frac{2\sqrt{2}T\eta}{n}\Big{)}.$

Similarly, the stability bound in Theorem 6 also implies

$\mathbb{E}_{S,\mathcal{A}}[F_{S}(\mathbf{w}^{*},\bar{\mathbf{v}}_{T})-F(\mathbf{w}^{*},\bar{\mathbf{v}}_{T})]\leq
4G^{2}(1+L/\rho)\Big{(}\frac{\sqrt{T}\eta}{\sqrt{n}}+\frac{2\sqrt{2}T\eta}{n}\Big{)}.$

According to Theorem C.3, we know

$\displaystyle\mathbb{E}_{S,\mathcal{A}}[R_{S}(\bar{\mathbf{w}}_{T})-F_{S}(\mathbf{w}^{*},\bar{\mathbf{v}}_{T})]$
$\displaystyle\leq\mathbb{E}_{\mathcal{A}}\big{[}\max_{\mathbf{v}\in\mathcal{V}}F_{S}(\bar{\mathbf{w}}_{T},\mathbf{v})-\min_{\mathbf{w}\in\mathcal{W}}F_{S}(\mathbf{w},\bar{\mathbf{v}}_{T})\big{]}$
$\displaystyle\leq
G^{2}\eta+\frac{D^{2}}{2T\eta}+\frac{2GK_{P}D+12G^{2}\eta\sum_{j=1}^{T}k_{j}+G\sum_{j=K_{P}}^{T}1/{j^{2}}}{T},$

where $D=D_{\mathbf{w}}+D_{\mathbf{v}}$. Putting the above three inequalities
back into Eq. (C.3), we obtain

$\displaystyle\mathbb{E}_{S,\mathcal{A}}[R(\bar{\mathbf{w}}_{T})-R(\mathbf{w}^{*})]\leq$
$\displaystyle
8G^{2}(1+L/\rho)\Big{(}\frac{\sqrt{T}\eta}{\sqrt{n}}+\frac{2\sqrt{2}T\eta}{n}\Big{)}+G^{2}\eta+\frac{D^{2}}{2T\eta}$
$\displaystyle+\frac{2GK_{P}D+12G^{2}\eta\sum_{j=1}^{T}k_{j}+G\sum_{j=K_{P}}^{T}1/{j^{2}}}{T}.$

If we choose $T\asymp n$ and $\eta\asymp(T\log(T))^{-1/2}$, combining the
above estimation with (C.15) implies

$\displaystyle\mathbb{E}_{S,\mathcal{A}}[R(\bar{\mathbf{w}}_{T})-R(\mathbf{w}^{*})]=\mathcal{O}\Big{(}\frac{(L/\rho)\log(n)}{\sqrt{n}\log(\lambda(P))}+\frac{K_{P}}{n\min\\{1,n\sqrt{\log(n)C_{P}\lambda(P)^{K_{P}}}\\}}\Big{)}.$

The above result combines with $K_{P}=0$ complete the proof. ∎

## References

  * [1] Shivani Agarwal and Partha Niyogi.  Generalization bounds for ranking algorithms via algorithmic stability.  Journal of Machine Learning Research, 10(Feb):441–474, 2009.
  * [2] Martin Arjovsky, Soumith Chintala, and Léon Bottou.  Wasserstein generative adversarial networks.  In International conference on machine learning, pages 214–223. PMLR, 2017.
  * [3] Yves F Atchadé, Gersende Fort, and Eric Moulines.  On perturbed proximal gradient algorithms.  The Journal of Machine Learning Research, 18(1):310–342, 2017.
  * [4] Ghadir Ayache and Salim El Rouayheb.  Private weighted random walk stochastic gradient descent.  IEEE Journal on Selected Areas in Information Theory, 2(1):452–463, 2021.
  * [5] Francis Bach and Eric Moulines.  Non-strongly-convex smooth stochastic approximation with convergence rate O(1/n).  In Advances in Neural Information Processing Systems, pages 773–781, 2013.
  * [6] Raef Bassily, Vitaly Feldman, Cristóbal Guzmán, and Kunal Talwar.  Stability of stochastic gradient descent on nonsmooth convex losses.  In Advances in Neural Information Processing Systems, 2020.
  * [7] Raef Bassily, Vitaly Feldman, Kunal Talwar, and Abhradeep Guha Thakurta.  Private stochastic convex optimization with optimal rates.  In Advances in Neural Information Processing Systems, pages 11279–11288, 2019.
  * [8] Raef Bassily, Kobbi Nissim, Adam Smith, Thomas Steinke, Uri Stemmer, and Jonathan Ullman.  Algorithmic stability for adaptive data analysis.  SIAM Journal on Computing, (50):377–405, 2021.
  * [9] Stéphane Boucheron, Gábor Lugosi, and Pascal Massart.  Concentration Inequalities: A Nonasymptotic Theory of Independence.  Oxford university press, 2013.
  * [10] Olivier Bousquet, Stéphane Boucheron, and Gábor Lugosi.  Introduction to statistical learning theory.  In Summer school on machine learning, pages 169–207. Springer, 2003\.
  * [11] Olivier Bousquet and André Elisseeff.  Stability and generalization.  The Journal of Machine Learning Research, 2:499–526, 2002.
  * [12] Olivier Bousquet, Yegor Klochkov, and Nikita Zhivotovskiy.  Sharper bounds for uniformly stable algorithms.  In Conference on Learning Theory, pages 610–626. PMLR, 2020.
  * [13] Zachary Charles and Dimitris Papailiopoulos.  Stability and generalization of learning algorithms that converge to global optima.  In International Conference on Machine Learning, pages 744–753, 2018.
  * [14] Jianshu Chen, Zaid J Towfic, and Ali H Sayed.  Dictionary learning over distributed models.  IEEE Transactions on Signal Processing, 63(4):1001–1016, 2014.
  * [15] Bo Dai, Albert Shaw, Lihong Li, Lin Xiao, Niao He, Zhen Liu, Jianshu Chen, and Le Song.  Sbeed: Convergent reinforcement learning with nonlinear function approximation.  In International Conference on Machine Learning, pages 1125–1134, 2018.
  * [16] Emily Diana, Wesley Gill, Michael Kearns, Krishnaram Kenthapadi, and Aaron Roth.  Minimax group fairness: Algorithms and experiments.  In Proceedings of the 2021 AAAI/ACM Conference on AI, Ethics, and Society, pages 66–76, 2021.
  * [17] Aymeric Dieuleveut and Francis Bach.  Nonparametric stochastic approximation with large step-sizes.  The Annals of Statistics, 44(4):1363–1399, 2016.
  * [18] Thinh T Doan, Lam M Nguyen, Nhan H Pham, and Justin Romberg.  Convergence rates of accelerated markov gradient descent with applications in reinforcement learning.  arXiv preprint arXiv:2002.02873, 2020.
  * [19] Thinh T Doan, Lam M Nguyen, Nhan H Pham, and Justin Romberg.  Finite-time analysis of stochastic gradient descent under markov randomness.  arXiv preprint arXiv:2003.10973, 2020.
  * [20] Simon S Du, Jianshu Chen, Lihong Li, Lin Xiao, and Dengyong Zhou.  Stochastic variance reduction methods for policy evaluation.  In International Conference on Machine Learning, pages 1049–1058, 2017.
  * [21] John Duchi and Yoram Singer.  Efficient online and batch learning using forward backward splitting.  In Advances in Neural Information Processing Systems, pages 495–503, 2009.
  * [22] John C Duchi, Alekh Agarwal, Mikael Johansson, and Michael I Jordan.  Ergodic mirror descent.  SIAM Journal on Optimization, 22(4):1549–1578, 2012.
  * [23] John C Duchi, Alekh Agarwal, and Martin J Wainwright.  Dual averaging for distributed optimization: Convergence analysis and network scaling.  IEEE Transactions on Automatic control, 57(3):592–606, 2011.
  * [24] Andre Elisseeff, Theodoros Evgeniou, and Massimiliano Pontil.  Stability of randomized learning algorithms.  Journal of Machine Learning Research, 6(Jan):55–79, 2005.
  * [25] Farzan Farnia and Asuman Ozdaglar.  Train simultaneously, generalize better: Stability of gradient-based minimax learners.  In International Conference on Machine Learning, pages 3174–3185. PMLR, 2021.
  * [26] Vitaly Feldman, Tomer Koren, and Kunal Talwar.  Private stochastic convex optimization: optimal rates in linear time.  In Proceedings of the 52nd Annual ACM SIGACT Symposium on Theory of Computing, pages 439–449, 2020.
  * [27] Vitaly Feldman and Jan Vondrak.  Generalization bounds for uniformly stable algorithms.  In Advances in Neural Information Processing Systems, pages 9747–9757, 2018.
  * [28] Vitaly Feldman and Jan Vondrak.  High probability generalization bounds for uniformly stable algorithms with nearly optimal rate.  In Conference on Learning Theory, pages 1270–1279. PMLR, 2019.
  * [29] Wei Gao, Rong Jin, Shenghuo Zhu, and Zhi-Hua Zhou.  One-pass auc optimization.  In International Conference on Machine Learning, pages 906–914, 2013.
  * [30] Ian Goodfellow, Jean Pouget-Abadie, Mehdi Mirza, Bing Xu, David Warde-Farley, Sherjil Ozair, Aaron Courville, and Yoshua Bengio.  Generative adversarial nets.  Advances in Neural Information Processing Systems, 27:2672–2680, 2014.
  * [31] Moritz Hardt and Tengyu Ma.  Identity matters in deep learning.  arXiv preprint arXiv:1611.04231, 2016.
  * [32] Moritz Hardt, Ben Recht, and Yoram Singer.  Train faster, generalize better: Stability of stochastic gradient descent.  In International Conference on Machine Learning, pages 1225–1234. PMLR, 2016.
  * [33] Bjorn Johansson, Maben Rabi, and Mikael Johansson.  A simple peer-to-peer algorithm for distributed optimization in sensor networks.  In 2007 46th IEEE Conference on Decision and Control, pages 4705–4710. IEEE, 2007.
  * [34] Björn Johansson, Maben Rabi, and Mikael Johansson.  A randomized incremental subgradient method for distributed optimization in networked systems.  SIAM Journal on Optimization, 20(3):1157–1170, 2010.
  * [35] Brian Kulis et al.  Metric learning: A survey.  Foundations and Trends® in Machine Learning, 5(4):287–364, 2013.
  * [36] Ilja Kuzborskij and Christoph Lampert.  Data-dependent stability of stochastic gradient descent.  In International Conference on Machine Learning, pages 2820–2829, 2018.
  * [37] Guanghui Lan.  First-order and Stochastic Optimization Methods for Machine Learning.  Springer Nature, 2020.
  * [38] Yunwen Lei, Zhenhuan Yang, Tianbao Yang, and Yiming Ying.  Stability and generalization of stochastic gradient methods for minimax problems.  In International Conference on Machine Learning, 2021.
  * [39] Yunwen Lei and Yiming Ying.  Fine-grained analysis of stability and generalization for stochastic gradient descent.  In International Conference on Machine Learning, pages 5809–5819. PMLR, 2020.
  * [40] Yunwen Lei and Yiming Ying.  Sharper generalization bounds for learning with gradient-dominated objective functions.  In International Conference on Learning Representations, 2020.
  * [41] Tian Li, Maziar Sanjabi, Ahmad Beirami, and Virginia Smith.  Fair resource allocation in federated learning.  arXiv preprint arXiv:1905.10497, 2019.
  * [42] Xiaoyu Li and Francesco Orabona.  On the convergence of stochastic gradient descent with adaptive stepsizes.  In The 22nd International Conference on Artificial Intelligence and Statistics, pages 983–992. PMLR, 2019.
  * [43] Junhong Lin and Lorenzo Rosasco.  Optimal learning for multi-pass stochastic gradient methods.  In Advances in Neural Information Processing Systems, pages 4556–4564, 2016.
  * [44] Tianyi Lin, Chi Jin, and Michael Jordan.  On gradient descent ascent for nonconvex-concave minimax problems.  In International Conference on Machine Learning, pages 6083–6093. PMLR, 2020.
  * [45] Mingrui Liu, Zhuoning Yuan, Yiming Ying, and Tianbao Yang.  Stochastic auc maximization with deep neural networks.  arXiv preprint arXiv:1908.10831, 2019.
  * [46] Mingrui Liu, Xiaoxuan Zhang, Zaiyi Chen, Xiaoyu Wang, and Tianbao Yang.  Fast stochastic AUC maximization with O(1/n)-convergence rate.  In International Conference on Machine Learning, pages 3195–3203, 2018.
  * [47] Tongliang Liu, Gábor Lugosi, Gergely Neu, and Dacheng Tao.  Algorithmic stability and hypothesis complexity.  In International Conference on Machine Learning, pages 2159–2167, 2017.
  * [48] Cassio G Lopes and Ali H Sayed.  Incremental adaptive strategies over distributed networks.  IEEE Transactions on Signal Processing, 55(8):4064–4077, 2007.
  * [49] Aleksander Madry, Aleksandar Makelov, Ludwig Schmidt, Dimitris Tsipras, and Adrian Vladu.  Towards deep learning models resistant to adversarial attacks.  arXiv preprint arXiv:1706.06083, 2017.
  * [50] Xianghui Mao, Kun Yuan, Yubin Hu, Yuantao Gu, Ali H Sayed, and Wotao Yin.  Walkman: A communication-efficient random-walk algorithm for decentralized optimization.  IEEE Transactions on Signal Processing, 68:2513–2528, 2020.
  * [51] Natalia Martinez, Martin Bertran, and Guillermo Sapiro.  Minimax pareto fairness: A multi objective perspective.  In International Conference on Machine Learning, pages 6755–6764. PMLR, 2020.
  * [52] Mehryar Mohri, Gary Sivek, and Ananda Theertha Suresh.  Agnostic federated learning.  In International Conference on Machine Learning, pages 4615–4625. PMLR, 2019.
  * [53] Angelia Nedic and Soomin Lee.  On stochastic subgradient mirror-descent algorithm with weighted averaging.  SIAM Journal on Optimization, 24(1):84–107, 2014.
  * [54] Arkadi Nemirovski, Anatoli Juditsky, Guanghui Lan, and Alexander Shapiro.  Robust stochastic approximation approach to stochastic programming.  SIAM Journal on Optimization, 19(4):1574–1609, 2009.
  * [55] Rufus Oldenburger.  Infinite powers of matrices and characteristic roots.  Duke Mathematical Journal, 6(2):357–361, 1940.
  * [56] Hassan Rafique, Mingrui Liu, Qihang Lin, and Tianbao Yang.  Weakly-convex–concave min–max optimization: provable algorithms and applications in machine learning.  Optimization Methods and Software, pages 1–35, 2021.
  * [57] Alexander Rakhlin, Ohad Shamir, and Karthik Sridharan.  Making gradient descent optimal for strongly convex stochastic optimization.  In Proceedings of the 29th International Coference on International Conference on Machine Learning, pages 1571–1578, 2012.
  * [58] S Sundhar Ram, A Nedić, and Venugopal V Veeravalli.  Incremental stochastic subgradient algorithms for convex optimization.  SIAM Journal on Optimization, 20(2):691–717, 2009.
  * [59] Herbert E Rauch, F Tung, and Charlotte T Striebel.  Maximum likelihood estimates of linear dynamic systems.  AIAA journal, 3(8):1445–1450, 1965.
  * [60] R Tyrrell Rockafellar.  Monotone operators and the proximal point algorithm.  SIAM Journal on Control and Optimization, 14(5):877–898, 1976.
  * [61] Ali H. Sayed.  Adaptive networks.  Proceedings of the IEEE, 102(4):460–497, 2014.
  * [62] Mark Schmidt, Nicolas Roux, and Francis Bach.  Convergence rates of inexact proximal-gradient methods for convex optimization.  In Advances in Neural Information Processing Systems, 2011.
  * [63] Suhail M Shah and Konstantin E Avrachenkov.  Linearly convergent asynchronous distributed admm via markov sampling.  arXiv preprint arXiv:1810.05067, 2018.
  * [64] Shai Shalev-Shwartz, Ohad Shamir, Nathan Srebro, and Karthik Sridharan.  Learnability, stability and uniform convergence.  Journal of Machine Learning Research, 11:2635–2670, 2010.
  * [65] Wei Shen, Zhenhuan Yang, Yiming Ying, and Xiaoming Yuan.  Stability and optimization error of stochastic gradient descent for pairwise learning.  Analysis and Applications, pages 1–41, 2019.
  * [66] Steve Smale and Ding-Xuan Zhou.  Online learning with markov sampling.  Analysis and Applications, 7(01):87–113, 2009.
  * [67] Padhraic Smyth.  Hidden markov models for fault detection in dynamic systems.  Pattern recognition, 27(1):149–164, 1994.
  * [68] Tao Sun and Dongsheng Li.  Decentralized markov chain gradient descent.  arXiv preprint arXiv:1909.10238, 2019.
  * [69] Tao Sun, Yuejiao Sun, and Wotao Yin.  On markov chain gradient descent.  In Advances in neural information processing systems, 2018.
  * [70] Vladislav B Tadić and Arnaud Doucet.  Asymptotic bias of stochastic gradient search.  The Annals of Applied Probability, 27(6):3255–3304, 2017.
  * [71] Enayat Ullah, Tung Mai, Anup Rao, Ryan Rossi, and Raman Arora.  Machine unlearning via algorithmic stability.  arXiv preprint arXiv:2102.13179, 2021.
  * [72] Vladimir Vapnik.  The nature of statistical learning theory.  Springer science & business media, 1999.
  * [73] Sharan Vaswani, Francis Bach, and Mark Schmidt.  Fast and faster convergence of sgd for over-parameterized models and an accelerated perceptron.  In The 22nd International Conference on Artificial Intelligence and Statistics, pages 1195–1204. PMLR, 2019.
  * [74] Puyu Wang, Yunwen Lei, Yiming Ying, and Hai Zhang.  Differentially private sgd with non-smooth losses.  Applied and Computational Harmonic Analysis, 56:306–336, 2022.
  * [75] Kilian Q Weinberger and Lawrence K Saul.  Distance metric learning for large margin nearest neighbor classification.  Journal of machine learning research, 10(2), 2009.
  * [76] Eric Xing, Michael Jordan, Stuart J Russell, and Andrew Ng.  Distance metric learning with application to clustering with side-information.  Advances in neural information processing systems, 15, 2002.
  * [77] Zhenhuan Yang, Yunwen Lei, Siwei Lyu, and Yiming Ying.  Stability and differential privacy of stochastic gradient descent for pairwise learning with non-smooth loss.  In International Conference on Artificial Intelligence and Statistics, pages 2026–2034. PMLR, 2021.
  * [78] Zhenhuan Yang, Yunwen Lei, Puyu Wang, Tianbao Yang, and Yiming Ying.  Simple stochastic and online gradient descent algorithms for pairwise learning.  Advances in Neural Information Processing Systems, 34, 2021.
  * [79] Yiming Ying and Peng Li.  Distance metric learning with eigenvalue optimization.  The Journal of Machine Learning Research, 13(1):1–26, 2012.
  * [80] Yiming Ying and Massimiliano Pontil.  Online gradient descent learning algorithms.  Foundations of Computational Mathematics, 8(5):561–596, 2008.
  * [81] Yiming Ying, Longyin Wen, and Siwei Lyu.  Stochastic online auc maximization.  Advances in neural information processing systems, 29, 2016.
  * [82] Yiming Ying and D-X Zhou.  Online regularized classification algorithms.  IEEE Transactions on Information Theory, 52(11):4775–4788, 2006\.
  * [83] Yiming Ying and Ding-Xuan Zhou.  Online pairwise learning algorithms.  Neural computation, 28(4):743–777, 2016.
  * [84] Sihan Zeng, Thinh T Doan, and Justin Romberg.  Finite-time analysis of decentralized stochastic approximation with applications in multi-agent and multi-task learning.  In 2021 60th IEEE Conference on Decision and Control (CDC), pages 2641–2646. IEEE, 2021.
  * [85] Junyu Zhang, Mingyi Hong, Mengdi Wang, and Shuzhong Zhang.  Generalization bounds for stochastic saddle point problems.  In International Conference on Artificial Intelligence and Statistics, pages 568–576. PMLR, 2021.
  * [86] Tong Zhang.  Solving large scale linear prediction problems using stochastic gradient descent algorithms.  In International Conference on Machine Learning, pages 919–926, 2004.
  * [87] Peilin Zhao, Steven CH Hoi, Rong Jin, and Tianbao Yang.  Online AUC maximization.  In International Conference on Machine Learning, pages 233–240, 2011.
  * [88] Yi Zhou, Yingbin Liang, and Huishuai Zhang.  Generalization error bounds with probabilistic guarantee for SGD in nonconvex optimization.  arXiv preprint arXiv:1802.06903, 2018.


# Non-thermal hot dark matter in light of the $S_{8}$ tension

Subinoy Das<EMAIL_ADDRESS>Indian Institute of Astrophysics, Bengaluru,
Karnataka 560034, India  Anshuman Maharana<EMAIL_ADDRESS>Harish-Chandra Research Institute, HBNI, Allahabad 211019, India.  Vivian
Poulin<EMAIL_ADDRESS>Laboratoire Univers & Particules de
Montpellier (LUPM), CNRS & Université de Montpellier (UMR-5299), Place Eugène
Bataillon, F-34095 Montpellier Cedex 05, France  Ravi Kumar Sharma
<EMAIL_ADDRESS>Indian Institute of Astrophysics, Bengaluru, Karnataka
560034, India

###### Abstract

The $\Lambda$CDM prediction of $S_{8}\equiv\sigma_{8}(\Omega_{m}/0.3)^{0.5}$ –
where $\sigma_{8}$ is the root mean square of matter fluctuations on a 8
$h^{-1}$Mpc scale – once calibrated on Planck CMB data is $2-3\sigma$ lower
than its direct estimate by a number of weak lensing surveys. In this paper,
we explore the possibility that the ’$S_{8}$-tension’ is due to a non-thermal
hot dark matter (HDM) fractional contribution to the universe energy density
leading to a power suppression at small-scales in the matter power spectrum.
Any HDM models can be characterized by its effective mass $m_{\rm sp}^{\rm
eff}$ and its contribution to the relativistic degrees of freedom at CMB
decoupling $\Delta N_{\rm eff}$. Taking the specific example of a sterile
particle produced from the decay of the inflaton during a matter dominated
era, we find that from Planck only the tension can be reduced below $2\sigma$,
but Planck does not favor a non-zero $\\{m_{\rm sp}^{\rm eff},\Delta N_{\rm
eff}\\}$. In combination with a measurement of $S_{8}$ from
KIDS1000+BOSS+2dfLenS, the $S_{8}$-tension would hint at the existence of a
particle of mass $m_{\rm sp}^{\rm eff}\simeq 0.67_{-0.48}^{+0.26}$ ${\rm eV}$
with a contribution to $\Delta N_{\rm eff}\simeq 0.06\pm 0.05$. However,
Pantheon and BOSS BAO/$f\sigma_{8}$ data restricts the particle mass to
$m_{\rm sp}^{\rm eff}\simeq 0.48_{-0.36}^{+0.17}$ and contribution to $\Delta
N_{\rm eff}\simeq 0.046_{-0.031}^{+0.004}$. We discuss implications of our
results for other canonical non-thermal HDM models– the Dodelson-Widrow model
and a thermal sterile particle with a different temperature in the hidden
sector. We report competitive results on such hidden sector temperature which
might have interesting implications for particle physics model building, in
particular connecting the $S_{8}$-tension to the longstanding short baseline
oscillation anomaly.

##  I Introduction

The $\Lambda$ Cold Dark Matter ($\Lambda$CDM) model of cosmology is incredibly
powerful at describing a wide variety of observation up to a high degree of
accuracy despite the nature of its dominant components - Cold Dark Matter
(CDM) and Dark Energy (DE) – still being unknown. Nevertheless, in recent
years, a number of intriguing discrepancies has emerged between the values of
some cosmological parameter predicted within $\Lambda$CDM – once the model is
calibrated onto Planck Cosmic Microwave Background (CMB) data, Baryon Acoustic
Oscillation (BAO) and luminosity distance to SuperNovae of type Ia (SNIa) –
and their direct measurements.

At the heart of this study is the longstanding tension affecting the
determination of the amplitude of the matter fluctuations, typically
parameterized as $S_{8}\equiv\sigma_{8}(\Omega_{m}/0.3)^{0.5}$, where
$\sigma_{8}$ is the root mean square of matter fluctuations on a 8 $h^{-1}$Mpc
scale, and $\Omega_{m}$ is the total matter abundance. The latest prediction
from Planck CMB data within the $\Lambda$CDM framework is $S_{8}=0.832\pm
0.013$ [1]. Originally, observations of galaxies through weak lensing by the
CFHTLenS collaboration have indicated that the $\Lambda$CDM model predicts a
$S_{8}$ value that is larger than the direct measurement at the $2\sigma$
level [2, 3]. This tension has since then been further established within the
KiDS/Viking data [4, 5], but is milder within the DES data [6]. However, a re-
analysis of the DES data, combined with KiDS/Viking, lead to a determination
of $S_{8}$ that is discrepant with Planck at the $3\sigma$ level,
$S_{8}=0.755_{-0.021}^{+0.019}$ [5]. Recently, the combination of KiDS/Viking
and SDSS data has established $S_{8}=0.766^{+0.02}_{-0.014}$ [7]. Moreover, it
is now understood that the tension is driven by a lower matter clustering
amplitude $\sigma_{8}$. This is particularly interesting for model building:
resolving the $S_{8}$ tension requires to decrease the amplitude of matter
fluctuations on scales $k\sim 0.1-1~{}h$/Mpc, which can be easily achieved in
a variety of models often related to new DM properties [8, 9, 10, 11, 12, 13,
14, 15, 16, 17, 18, 19, 20], or new neutrino properties [21, 22].

In this paper, we explore the possibility that the ‘$S_{8}$-tension’ is due to
the existence of a non-thermal hot dark matter (HDM) particle contributing to
a fraction of the dark matter (DM) density in the universe and leading to a
power suppression at small-scales in the matter power spectrum. It is well
known that just adding a thermal neutrino like radiation $\Delta N_{\rm eff}$,
together with a non-zero neutrino mass $m_{\nu}$, does not resolve the $S_{8}$
problem [21, 1]. Here, we explore the consequences of a non-thermal momentum
distribution for the hot component (or a temperature different from our
visible sector), for the $S_{8}$-tension. In practice, we consider the
momentum distribution associated with sterile particle produced from decays
during a matter domination to radiation domination transition of the universe.
We dub the model $\nu_{\rm NT}\Lambda$CDM. From the point of view of
theoretical models, it is natural for the early universe to enter an epoch of
matter to radiation transition (this happens during the decay of the inflation
or if cold moduli dominates the energy density of the universe at early
times). Such particles have a characteristic momentum distribution, that is
associated with decays taking place in a matter dominated universe evolving to
radiation domination.

For decay products that thermalize, thermalization leads to loss of all
information about the kinematics of the decay process. But in a setting with a
large number of hidden sectors, one can expect that some of the species
produced during the decay do not thermalize due to very weak interactions. Our
scenario belongs to such category where a moduli or inlfaton field decays to
non-thermal sterile particle. There might be other particle like Feebly
interacting massive particle (FIMP) which can also produce non-thermal or
partially thermal neutrino like particle [23]. The presence of a non-thermal
dark radiation can affect CMB [24] as well as large scale structure in
specific ways and can be probed by precision cosmological data. The study of
the implications of sterile particles with this momentum distribution for
precision cosmology was recently initiated111For earlier in work on inert
particles from decays see e.g [25, 26, 27] in Ref. [28]. Given that the effect
of massive sterile particles on the CMB and matter power spectra is well-known
(e.g. [29, 30, 31] for reviews), it was anticipated there one might get a
strong power suppression in the matter power spectra due to the typical
momentum distribution of non-thermal decay product. This power suppression has
implications for the $S_{8}$-tension.

In this article, we perform a comprehensive monte-carlo markov chains (MCMC)
analysis against up-to-date data from Planck, BOSS (BAO and $f\sigma_{8}$) and
Pantheon data, with and without the inclusion of a prior on the value of
$S_{8}$ as measured with the KiDS/Viking+BOSS+2dFLens data222For analysis in
similar spirit (although without inclusion of the $S_{8}$) prior motivated by
short base line neutrino experiments see e.g [32, 33]. Of course, here the
momentum distribution of the sterile particles is assumed to be as motivated
by neutrino physics i.e thermal or the Dodelson-Widrow distribution [34].. We
find that the $\nu_{\rm NT}\Lambda$CDM model can indeed alleviate the tension
between Planck and $S_{8}$ measurements, but the success of the resolution is
degraded once BOSS and Pantheon data are included in the analysis. To better
understand the model features leading to a resolution of the tension, we
compare the non-thermal sterile neutrino model to the standard massive
neutrino model with extra relativistic degrees of freedom. We find that, for
similar effect on the CMB power spectrum, the $\nu_{\rm NT}\Lambda$CDM leads
to a much stronger suppression in the matter power spectrum at late-times, and
therefore to a more significant decrease in $\sigma_{8}$. The impact of the
$\nu_{\rm NT}\Lambda$CDM is barely visible on the BAO scale and luminosity
distance, but does significantly affect $f\sigma_{8}$ predictions. The model
is therefore further constrained by BOSS data. Future measurements of the
matter power-spectrum and $f\sigma_{8}$ at late-times will further test this
scenario.

Although the MCMC analysis is carried out for sterile particles with the above
described momentum distributions, it has implications for a wide class of
models. As is well known (see e.g. [29, 35]), the cosmological implications of
a hot a sterile component is captured effectively by just two parameters: a)
the contribution of the component to the present day energy density, usually
reported in terms of the effective mass parameter: $m_{\rm eff}$ b) the
contribution of the component to the energy density at the time of the CMB
decoupling, usually reported in terms of $\Delta N_{\rm eff}$ 333In the models
we will discussing, this is same as $\Delta N_{\rm eff}$ at the time of
neutrino decoupling . These parameters are determined by the first two moments
of the momentum distribution and the mass of the sterile particle. Two models
with equal values of $m_{\rm eff}$ and $\Delta N_{\rm eff}$ will have the same
physics even if the functional form of the associated momentum distributions
is different. We translate the results of the analysis for our model
parameters to results on the effective parameters. Our results therefore have
direct implications for other well motivated momentum distributions such as a
thermal distribution with a different temperature from that the of Standard
Model[36, 37, 38], the Dodelson-Widrow distribution [34] or distributions
similar to the Dodeldon-Widrow discussed in Refs .[39, 40].

Our paper is structured as follows: In Sec. II, we present our model and the
mapping onto generic phenomenological parameters; in Sec. III, we perform a
MCMC analysis against a suite of up-to-date cosmological data and discuss the
extent at which the $\nu_{\rm NT}\Lambda$CDM can resolve the $S_{8}$-tension;
in Sec. IV, we draw implications of our results for other HDM models; finally,
we conclude in Sec. V.

##  II Non-thermal hot dark matter

###  II.1 The model

The physics of a constituent species of dark matter depends on its mass,
interactions and also on its momentum distribution function. For species that
thermalise, the process of thermalisation brings the momentum distribution to
the Fermi-Dirac or Bose-Einstein form. On the other hand, for non-thermal
constituents the momentum distribution is determined by their production
mechanism. Thus, it is important to isolate natural production mechanisms for
species that can constitute dark matter, the associated momentum distributions
and their implications for cosmology.

In this section, we will review the basics of the production mechanism and the
form of the momentum distribution that we will be considering. Our discussion
will be brief, we refer the reader to [28] and the references therein for the
details. At early times, the energy density of the universe will be taken to
be dominated by cold particles of a species $\varphi$. We will denote the mass
of the particles of $\varphi$ by $m_{\varphi}$ and their decay width to be
$\tau$. We will be focusing on the case when the $\varphi$ is the inflaton,
with inflation taking place at the GUT scale and decays of the inflaton taking
place due to a non-renormalizable interaction at the GUT scale. Thus we will
take $m_{\varphi}\sim 10^{-6}M_{\rm pl}$ and $\tau\sim 10^{8}/m_{\varphi}$.
The branching ratio of the $\varphi$ particles to the sterile particles will
be taken to be $B_{\rm sp}$, the sterile particles so produced will be taken
not to thermalise. We will assume that the other decay products thermalise, as
this sector would contain the Standard Model, we will refer to it as the
Standard Model sector. All decay products will be taken to be relativistic at
the time of production. As the $\varphi$ particles decay, the universe goes
into a matter to radiation epoch, finally becoming fully radiation dominated.

During the matter to radiation dominated epoch the evolution of the universe
is governed by the equations

$\dot{\rho}_{\rm mat}+3H\rho_{\rm mat}=-{\rho_{\rm mat}\over\tau},$ (1)
$\dot{\rho}_{\rm rad}+4H\rho_{\rm rad}=+{\rho_{\rm mat}\over\tau},$ (2)

and

$H=\left({\dot{a}\over a}\right)=\sqrt{{{\rho_{\rm mat}+\rho_{\rm rad}}\over
3M^{2}_{\rm pl}}}.$ (3)

In the above, $\rho_{\rm mat}$ denotes the energy density in matter and
$\rho_{\rm rad}$ is the energy density in radiation. The energy density in
radiation is the sum of the energy densities in the Standard Model sector and
the sterile particles (since the sterile particles are highly relativistic at
the time of production, thus they contribute to the energy density as
radiation during the epoch that the decays take place). It is useful to
introduce the dimensionless variables

$\displaystyle\theta={t\over\tau},\ \ \ \hat{s}(\theta)=a(\tau\theta),\ \ \ $
$\displaystyle e_{\rm mat}(\theta)={\tau^{2}\rho_{\rm mat}(\tau\theta)\over
M_{\rm pl}^{2}}\ \ \ \textrm{and}\ \ \ e_{\rm rad}(\theta)={\tau^{2}\rho_{\rm
rad}(\tau\theta)\over M_{\rm pl}^{2}}.\ \ \ $ (4)

Once almost all $\varphi$ particles have decayed, one can take the universe to
be composed of a thermal bath (which contains the Standard Model sector) and
the sterile particles governed by the standard cosmological evolution
equations. In practice, we will start with a matter dominated universe at an
‘initial time’ $(t=\theta=0)$, evolve the universe using the equations 1, 2
and (3) up to a fiducial dimensionless time $\theta^{*}$ which is large enough
so that a good fraction the $\varphi$ particles have decayed by that time (we
will choose $\theta^{*}=15$ for this) and use the results of this as the
initial conditions for the standard cosmological evolution.

The momentum distribution of the sterile particles can be computed using the
following: as a result of the decays the co-moving number density of the
sterile particles falls off as $N(t)=N(0)e^{-t/\tau}$ (and the branching ratio
to the sterile particles is $B_{\rm sp}$), once produced the sterile particles
free-stream. We will be making use of publicly available package CLASS [41,
42, 43] to incorporate the effects of the sterile particles, the input
required for this the momentum distribution of the sterile particles today.
This was obtained in [28] to be

$f(q)={32\over{\pi\hat{E}^{3}}}\left({N(0)B_{\rm
sp}\over{\hat{s}^{3}(\theta^{*})}}\right){e^{-\hat{s}^{-1}(y)}\over
q^{3}\hat{H}(\hat{s}^{-1}(y))},$ (5)

where

$y={q\over 4}\hat{s}(\theta^{*}),$ (6)

and the range of q is given by

${4\over\hat{s}(\theta^{*})}<q<4.$ (7)

where $\hat{E}=m_{\varphi}/2$, $N(0)$ is the initial number density of the
$\varphi$ particles, $\hat{s}(\theta^{*})$ the value scale factor444In (5) the
convention is that the scale factor is 1 at $\theta=0$ at the fiducial
dimensionless time $\theta^{*}$, $\hat{s}^{-1}$ the functional inverse of the
scale factor function as a function of the dimensionless time and
$\hat{H}=\hat{s}^{\prime}(\theta)/\hat{s}(\theta)$ the dimensionless Hubble
constant. The momentum $q$ in (5) is the magnitude of the momentum in units of
the typical momentum of the sterile particles today $(T_{\rm ncdm,0})$. The
typical momentum was found to be

$T_{\rm ncdm,0}=0.418\left(m_{\varphi}^{2}\tau\over M_{\rm
pl}\right)^{1/2}{T_{\rm cmb}\over{(1-B_{\rm sp})^{1/4}}}\equiv\zeta T_{\rm
cmb}$ (8)

in [28]. The distribution function in (5) is in units of $T_{\rm ncdm,0}^{3}$.
Note that although naively $f(q)$ seems to depend on $N(0)$, the full
expression is independent of $N(0)$ as long as we take the universe to be
completely matter dominated at the initial times. It is interesting to compare
the distribution to a thermal one, as shown in figure 1. For the same value of
$\Delta N_{\rm eff}$, the distribution is peaked at higher values of the
momentum, but is much broader. The mean momentum of sterile particles is
greater than that of the CMB by the factor $\zeta$ defined in (8). For our
choice of parameters $\zeta\sim 5$. The sterile particles become non-
relativistic when their typical momentum becomes of the order of their mass
i.e the temperature of the Standard Model plasma becomes of the order
${m_{\varphi}/5}$.

Figure 1: Comparative plot with a thermal distribution. The nonthermal
distribution is plotted in orange and is for value of the parameter
$m_{\varphi}=10^{-6}M_{\rm pl}$, $\tau={10^{8}/m_{\varphi}}$ . The thermal
distribution is in blue. The momenta and the distribution functions for both
plots are in units of $T_{\rm ncdm,0}$ as associated with the above values of
$m_{\varphi}$ and $\tau$. $\Delta N_{\rm eff}$ is taken to be 0.15 for both
distributions.

###  II.2 Mapping onto generic parameters

Our model has four microscopic parameters: $m_{\varphi},\tau$ (the mass and
lifetime of the decaying particle), $B_{\rm sp}$ (the branching ratio for
decay to the sterile particle) and $m_{\rm sp}$ (the mass of the sterile
particle) in addition to those of $\Lambda$CDM. Our choice of the first two
parameters ($m_{\varphi}\sim 10^{-6}M_{\rm pl}$ and $\tau\sim
10^{8}/m_{\varphi}$) is motivated by taking $\varphi$ to be driving inflation
at the GUT scale and decaying by GUT scale interactions. On the other hand,
the other parameters $B_{\rm sp}$ and $m_{\rm sp}$ will be traded for
effective parameters more directly connected to observables. Indeed, the
physical effects of new sterile particle/species on the cosmological
background and perturbation evolution can be completely described by three
parameters: $\Delta N_{\rm eff}$ (the effective number of relativistic
neutrinos at the time of neutrino decoupling), $w_{\rm sp}\equiv\Omega_{\rm
sp}h^{2}$ (where $\Omega_{\rm sp}$ is the fractional contribution of the
particle to today’s energy density and h the reduced Hubble parameter, this is
often characterised by the effective mass of the particle $m_{\rm sp}^{\rm
eff}.=w_{\rm sp}94.05{\rm eV}$) and $\lambda_{\rm FS}$ (the free-streaming
length associated with species). The free-streaming length is determined once
the first two quantities are known, hence effectively there are two parameters
[29]. Here we can take $\Delta N_{\rm eff}$ and $m_{\rm sp}^{\rm eff}$ as two
independent quantities, defined as

$\Delta N_{\rm eff}\equiv\frac{\rho_{s}^{\rm
rel}}{\rho_{\nu}}=\left[{\frac{1}{\pi^{2}}\int\\!\\!dp\,\,p^{3}f(p)}\right]/\left[{\frac{7}{8}\frac{\pi^{2}}{15}{T_{\nu}^{\rm
id}}^{4}}\right]$ (9)

with $T_{\nu}^{\rm id}\equiv(4/11)^{1/3}T_{\gamma}$ and

$\frac{m_{\rm sp}^{\rm eff}}{94.05{\rm
eV}}\equiv\omega_{s}\equiv\Omega_{s}h^{2}=\left[{\frac{m_{\rm
sp}}{\pi^{2}}\int\\!\\!dp\,\,p^{2}f(p)}\right]\times\left[\frac{h^{2}}{\rho_{c}^{0}}\right]$
(10)

where $\rho_{c}^{0}$ is the critical density today and $h$ the reduced Hubble
parameter. In our model, the effective parameters $m_{\rm eff}$ and $\Delta
N_{\rm eff}$ in terms of the microscopic parameters are given by [28]:

$\Delta N_{\rm eff}={43\over 7}{B_{\rm sp}\over{1-B_{\rm
sp}}}\left({g_{*}(T(t_{\nu}))\over g_{*}(T(t^{*}))}\right)^{1/3}$ (11)

and

$m^{\rm eff}_{\rm sp}={62.1m_{\rm sp}\over{g_{*}^{1/4}(T(t^{*}))}}{B_{\rm
sp}\over{(1-B_{\rm sp})^{3/4}}}\left({M_{\rm pl}\over\tau
m_{\varphi}^{2}}\right)^{1/2},$ (12)

where $g_{*}(T(t_{\nu}))$ and $g_{*}(T(t^{*}))$ are the effective number of
degrees of freedom at the time of neutrino decoupling and the end of the
reheating epoch (we will take the later to be equal to $100$). We will thus
scan over $m_{\rm sp}$ and $B_{\rm sp}$, and use Eqs. 12 and 11 to relate to
phenomenological parameters. In sec. IV, we will then translate our results
for two other models of interest.

Model |  $\Lambda$CDM |  $\nu\Lambda$CDM |  $\nu_{\rm NT}\Lambda$CDM
---|---|---|---
Parameter | Planck | Planck + $S_{8}$ | Planck | Planck + $S_{8}$ | Planck | Planck + $S_{8}$
$100~{}\omega_{b}$ | $2.24(2.24)_{-0.015}^{+0.014}$ | $2.252(2.256)_{-0.015}^{+0.014}$ | $2.244_{-0.018}^{+0.016}$ | $2.257(2.259)\pm 0.017$ | $2.242(2.247)_{-0.016}^{+0.014}$ | $2.247(2.247)_{-0.015}^{+0.014}$
$\omega_{\rm cdm}$ | $0.1198(0.1195)_{-0.0012}^{+0.0013}$ | $0.1182(0.1177)\pm 0.0011$ | $0.1217_{-0.002}^{+0.0015}$ | $0.1198(0.1182)_{-0.0018}^{+0.0013}$ | $0.1182(0.1198)_{-0.0021}^{+0.0039}$ | $0.1142(0.1110)_{-0.003}^{+0.0049}$
$100*\theta_{s}$ | $1.04190(1.04178)_{-0.0003}^{+0.00029}$ | $1.04202(1.04217)_{-0.0003}^{+0.00029}$ | $1.04166_{-0.00033}^{+0.00037}$ | $1.04179(1.04191)_{-0.00032}^{+0.00035}$ | $1.04180(1.04187)\pm 0.00031$ | $1.04186(1.04190)_{-0.00029}^{+0.00031}$
$n_{s}$ | $0.9661(0.9663)_{-0.0043}^{+0.0041}$ | $0.9695(0.971)_{-0.0041}^{+0.0039}$ | $0.9685_{-0.006}^{+0.0049}$ | $0.9717(0.9732)_{-0.0056}^{+0.0048}$ | $0.9652(0.9677)_{-0.0051}^{+0.0044}$ | $0.9652(0.9661)_{-0.0045}^{+0.0047}$
${\rm ln}(10^{10}A_{s})$ | $3.044(3.044)\pm 0.014$ | $3.041(3.042)_{-0.015}^{+0.014}$ | $3.052_{-0.016}^{+0.015}$ | $3.048(3.050)_{-0.017}^{+0.016}$ | $3.048(3.0480)\pm 0.015$ | $3.046(3.044)_{-0.016}^{+0.014}$
$\tau_{\rm reio}$ | $0.0541(0.0541)_{-0.0071}^{+0.0075}$ | $0.0542(0.0556)_{-0.0078}^{+0.0074}$ | $0.0558_{-0.0081}^{+0.0073}$ | $0.0555(0.0590)_{-0.0082}^{+0.0077}$ | $0.0546(0.0559)_{-0.0078}^{+0.0072}$ | $0.0548(0.0536)_{-0.0079}^{+0.0069}$
$m_{\nu}$ [eV] | $-$ | $-$ | $<0.073$ | $<0.1(0)$ | $-$ | $-$
$m_{\rm sp}^{\rm eff}$ [eV] | $-$ | $-$ | $-$ | $-$ | $<0.93(0)$ | $0.67(0.90)_{-0.48}^{+0.26}$
$\Delta N_{\rm eff}$ | $-$ | $-$ | $<0.28$ | $<0.24(0.03)$ | $<0.18(0.03)$ | $0.0614(0.034)_{-0.047}^{+0.0052}$
$S_{8}$ | $0.834(0.832)\pm 0.013$ | $0.814(0.809)_{-0.011}^{+0.01}$ | $0.834(0.838)_{-0.013}^{+0.013}$ | $0.812(0.814)\pm 0.011$ | $0.816(0.831)_{-0.016}^{+0.022}$ | $0.789(0.791)\pm 0.016$
$\Omega_{m}$ | $0.3078(0.3068)_{-0.0076}^{+0.0074}$ | $0.2981(0.2948)_{-0.0066}^{+0.0061}$ | $0.3154(0.3084)_{-0.015}^{+0.0094}$ | $0.3084(0.295)_{-0.018}^{+0.0081}$ | $0.314(0.305)_{-0.01}^{+0.008}$ | $0.311(0.308)_{-0.01}^{+0.008}$
$H_{0}$ [km/s/Mpc] | $68(68.04)\pm 0.56$ | $68.73(68.99)_{-0.51}^{+0.49}$ | $67.83(67.95)_{-1}^{+1.2}$ | $68.26(69.11)_{-0.93}^{+1.5}$ | $67.75(68.34)_{-0.65}^{+0.62}$ | $67.91(68.04)_{-0.61}^{+0.67}$
$\chi^{2}_{\rm min}$ | 2774.8 | 2783.4 | 2774.9 | 2782.0 | 2775.0 | 2778.60

Table 1:  The mean (best-fit) $\pm 1\sigma$ error of the cosmological
parameters in the $\Lambda$CDM and $\nu_{\rm NT}\Lambda$CDM model obtained
from the analysis of Planck [44] and Planck+$S_{8}$ [7] data. The definition
of $m_{\rm sp}^{\rm eff}$ is given in Eq. 12. Upper limits are given at the
95% C.L.

Model |  $\Lambda$CDM |  $\nu\Lambda$CDM |  $\nu_{\rm NT}\Lambda$CDM
---|---|---|---
Parameter | Planck+Ext | Planck+Ext+$S_{8}$ | Planck+Ext | Planck+Ext+$S_{8}$ | Planck+Ext | Planck+Ext+$S_{8}$
$100~{}\omega_{b}$ | $2.241(2.238)_{-0.014}^{+0.013}$ | $2.248(2.258)\pm{0.013}$ | $2.249(2.248)\pm 0.015$ | $2.257(2.250)\pm 0.015$ | $2.245(2.245)\pm 0.014$ | $2.250(2.253)_{-0.014}^{+0.013}$
$\omega_{\rm cdm}$ | $0.1197(0.1204)\pm 0.0009$ | $0.1187(0.1182)_{-0.0008}^{+0.0009}$ | $0.121(0.1194)_{-0.0019}^{+0.0012}$ | $0.1198(0.1186)_{-0.0017}^{+0.0011}$ | $0.1179(0.1179)_{-0.0018}^{+0.0031}$ | $0.1152(0.1101)_{-0.0023}^{+0.0036}$
$100*\theta_{s}$ | $1.04192(1.04204)_{-0.00029}^{+0.00028}$ | $1.04197(1.04186)_{-0.00029}^{+0.0003}$ | $1.04172(1.04194)_{-0.00031}^{+0.00034}$ | $1.04179(1.04194)_{-0.00031}^{+0.00036}$ | $1.04188(1.04193)_{-0.00029}^{+0.0003}$ | $1.04193(1.04194)_{-0.00028}^{+0.00029}$
$n_{s}$ | $0.9664(0.9660)_{-0.0037}^{+0.0038}$ | $0.9683(0.9705)_{-0.0038}^{+0.0036}$ | $0.9699(0.9693)_{-0.0049}^{+0.0044}$ | $0.9721(0.9706)_{-0.0048}^{+0.0043}$ | $0.9666(0.9664)_{-0.0041}^{+0.0039}$ | $0.9669(0.9678)_{-0.004}^{+0.0039}$
${\rm ln}(10^{10}A_{s})$ | $3.044(3.05)_{-0.015}^{+0.014}$ | $3.038(3.045)_{-0.015}^{+0.013}$ | $3.052(3.049)_{-0.016}^{+0.014}$ | $3.046(3.035)_{-0.016}^{+0.015}$ | $3.048(3.052)_{-0.015}^{+0.014}$ | $3.046(3.054)_{-0.015}^{+0.014}$
$\tau_{\rm reio}$ | $0.0542(0.0574)_{-0.0073}^{+0.0069}$ | $0.0526(0.056)_{-0.0076}^{+0.0069}$ | $0.0561(0.0569)_{-0.0081}^{+0.0066}$ | $0.0548(0.0515)_{-0.0081}^{+0.0073}$ | $0.0558(0.0576)_{-0.0076}^{+0.007}$ | $0.0556(0.0586)_{-0.0076}^{+0.0068}$
$m_{\nu}$ [eV] | $-$ | $-$ | $<0.040(0.005)$ | $<0.057(0.01)$ | $-$ | $-$
$m_{\rm sp}^{\rm eff}$ [eV] | $-$ | $-$ | $-$ | $-$ | $<0.73(0.21)$ | $0.48(0.92)_{-0.36}^{+0.17}$
$\Delta N_{\rm eff}$ | $-$ | $-$ | $<0.27(0.02)$ | $<0.26(0.006)$ | $<0.12(0.02)$ | $0.0457(0.0336)_{-0.031}^{+0.0038}$
$S_{8}$ | $0.832(0.842)\pm 0.011$ | $0.818(0.815)_{-0.0094}^{+0.0091}$ | $0.830(0.827)\pm 0.011$ | $0.814(0.815)_{-0.0097}^{+0.01}$ | $0.814(0.820)_{-0.014}^{+0.017}$ | $0.795(0.787)_{-0.013}^{+0.015}$
$\Omega_{m}$ | $0.3067(0.31)\pm 0.0055$ | $0.3007(0.2974)_{-0.0049}^{+0.0051}$ | $0.3084(0.3042)_{-0.006}^{+0.0059}$ | $0.3045(0.3037)_{-0.0072}^{+0.0061}$ | $0.309(0.308)_{-0.0059}^{+0.0059}$ | $0.306(0.304)\pm 0.006$
$H_{0}$ [km/s/Mpc] | $68.07(67.82)_{-0.43}^{+0.41}$ | $68.52(68.78)_{-0.4}^{+0.38}$ | $68.35(68.33)_{-0.7}^{+0.56}$ | $68.58(68.28)_{-0.73}^{+0.64}$ | $68.06(67.97)_{-0.48}^{+0.42}$ | $68.22(68.37)_{-0.43}^{+0.41}$
$\chi^{2}_{\rm min}$ | 3810.4 | 3818.2 | 3809.5 | 3816.4 | 3809.7 | 3814.5

Table 2: Same as Tab.1 with the addition of ’Ext’ data, which refers to the
combination BAO/FS+Pantheon.

##  III Resolving the $S_{8}$ tension with a non-thermal sterile neutrino

###  III.1 Details of the analysis

We perform a comprehensive MonteCarlo Markov Chain (MCMC) analysis and
confront the non thermal hot dark matter model to various combination of the
following data-sets:

  * •

Planck 2018 measurements of the low-$\ell$ CMB TT, EE, and high-$\ell$ TT, TE,
EE power spectra, together with the gravitational lensing potential
reconstruction [44].

  * •

The BAO measurements from 6dFGS at $z=0.106$ [45], SDSS DR7 at $z=0.15$ [46],
BOSS DR12 at $z=0.38,0.51$ and $0.61$ [47], and the joint constraints from
eBOSS DR14 Ly-$\alpha$ auto-correlation at $z=2.34$ [48] and cross-correlation
at $z=2.35$ [49].

  * •

The measurements of the growth function $f\sigma_{8}(z)$ (FS) from the CMASS
and LOWZ galaxy samples of BOSS DR12 at $z=0.38$, $0.51$, and $0.61$ [47].

  * •

The Pantheon SNIa catalogue, spanning redshifts $0.01<z<2.3$ [50].

  * •

The KIDS1000+BOSS+2dfLenS weak lensing data, compressed as a a split-normal
likelihood on the parameter $S_{8}=0.766^{+0.02}_{-0.014}$ [7].

Our baseline cosmology consists in the following combination of the six
$\Lambda$CDM parameters $\\{\omega_{b},\omega_{\rm
cdm},100\times\theta_{s},n_{s},{\rm ln}(10^{10}A_{s}),\tau_{\rm reio}\\}$,
plus two parameters describing the non-thermal hot dark matter, namely
$\\{m_{\rm sp},B_{\rm sp}\\}$. We dub this model $\nu_{\rm NT}\Lambda$CDM.
Standard model neutrinos are assumed to be massless.

To better understand how the $\nu_{\rm NT}\Lambda$CDM model can resolve the
$S_{8}$ tension, we will compare it to the standard $\Lambda$CDM model with
massless neutrinos, as well as to the $\Lambda$CDM model with free neutrino
masses $m_{\nu}$ and additional relativistic degrees of freedom $\Delta N_{\rm
eff}$. In that latter case, we assume degenerate neutrino masses and a free-
streaming $\Delta N_{\rm eff}$. We dub this model $\nu\Lambda$CDM. We run our
MCMCs with the Metropolis-Hasting algorithm as implemented in the
MontePython-v3 [51] code interfaced with our modified version of CLASS. All
reported $\chi^{2}_{\rm min}$ are obtained with the python package iMinuit
555https://iminuit.readthedocs.io/ [52]. We make use of a Choleski
decomposition to better handle the large number of nuisance parameters [53]
and consider chains to be converged with the Gelman-Rubin convergence
criterium $R-1\lesssim 0.05$ [54].

Figure 2: Reconstructed 2D posterior distributions of $\\{m_{\rm sp}^{\rm
eff},\Delta N_{\rm eff},S_{8},\Omega_{m}\\}$ with Planck and Planck+$S_{8}$
data (left panel) or Planck+BAO+SN1a and Planck+BAO+SN1a+$S_{8}$ data (right
panel).

Figure 3: Same as fig. 2 in the thermal neutrino case.

###  III.2 Results

We run two sets of analysis; in the first one, we confront the $\Lambda$CDM,
$\nu\Lambda$CDM and $\nu_{\rm NT}\Lambda$CDM models to Planck only and
Planck+$S_{8}$. In the second one, we add the BAO and Pantheon data to our
analysis. Our main results are reported in Tabs. 1 and 2 and displayed on
Figs. 2 and 3. We report results in the $\nu_{\rm NT}\Lambda$CDM in terms of
$\Delta N_{\rm eff}$ and $m_{\rm sp}^{\rm eff}$ defined in Eqs. 11 and 12. We
give the $\chi^{2}_{\rm min}$ per experiment in App. A.

####  III.2.1 Planck only

When the $\nu_{\rm NT}\Lambda$CDM model is confronted to Planck only, we
obtain a bound666Hereinafter, we quote 1-sided constraints at 95%C.L., and
two-sided ones at 68%C.L. on the mass $m_{\rm sp}^{\rm eff}<0.93$ eV and
$\Delta N_{\rm eff}<0.18$. Similarly, in the $\nu\Lambda$CDM case we obtain
$m_{\nu}<0.073$ eV and $\Delta N_{\rm eff}<0.28$ (recall that this limit
applies to individual neutrino masses in the degenerate case). The
$\chi^{2}_{\rm min}$ of Planck in the $\nu\Lambda$CDM and $\nu_{\rm
NT}\Lambda$CDM scenario is not improved over that of $\Lambda$CDM. We note
that the $\nu_{\rm NT}\Lambda$CDM model predicts a lower $S_{8}$ value than
other models. Indeed, we find $S_{8}(\nu\Lambda{\rm
CDM})=0.831_{-0.013}^{+0.012}$ and $S_{8}(\Lambda{\rm CDM})=0.832\pm 0.011$,
to be compared to $S_{8}(\nu_{\rm NT}\Lambda{\rm
CDM})=0.816_{-0.016}^{+0.022}$, i.e., a $\gtrsim 1\sigma$ downward shift. As a
result, the $S_{8}$ tension is alleviated from the $\sim 2.7\sigma$ level to
the $\sim 1.9\sigma$ level in the non-thermal HDM model. We note that our
constraints on $\Delta N_{\rm eff}$ in the non-thermal case is stronger than
that reported in Ref. [44] (constraints are identical in the thermal case).
This likely comes from the impact of running on physical parameters as opposed
to phenomenological parameters when exploring the parameter space.

Including the prior on $S_{8}$, we notice a mild detection of non-zero $m_{\rm
sp}^{\rm eff}=0.67_{-0.48}^{+0.26}$ eV and $\Delta N_{\rm
eff}=0.0614_{-0.047}^{+0.0052}$ in the $\nu_{\rm NT}\Lambda{\rm CDM}$ model,
while the constraints on the thermal neutrino mass simply relaxes to
$m_{\nu}<0.1$ eV. This translates into a reconstructed $S_{8}(\nu_{\rm
NT}\Lambda{\rm CDM})=0.789\pm 0.016$ and $S_{8}(\nu\Lambda{\rm CDM})=0.812\pm
0.011$, to be compared with the baseline $S_{8}(\Lambda{\rm
CDM})=0.814_{-0.011}^{+0.01}$. As a consequence, the $\chi^{2}_{\rm min}$ in
the combined analysis is lower in the non-thermal HDM case
$\Delta\chi^{2}_{\rm min}(\nu_{\rm NT}\Lambda{\rm CDM})=\chi^{2}_{\rm
min}(\Lambda{\rm CDM})-\chi^{2}_{\rm min}(\nu_{\rm NT}\Lambda{\rm CDM})=-4.8$
than in the thermal neutrino case $\Delta\chi^{2}_{\rm min}(\nu\Lambda{\rm
CDM})=\chi^{2}_{\rm min}(\Lambda{\rm CDM})-\chi^{2}_{\rm min}(\nu\Lambda{\rm
CDM})=-1.4$. If the $S_{8}$ tension worsens in the future, it would be
interesting to perform a more complete Bayesian analysis comparing these
models. We notice however that the total $\chi^{2}_{\rm min}$ is much less
significantly affected by the inclusion of the $S_{8}$ prior in the non-
thermal case (+3.6) than in the thermal case (+7.9), which is encouraging and
indicates that the $\nu_{\rm NT}\Lambda{\rm CDM}$ model can potentially
alleviate the tension between Planck and KIDS+BOSS. It remains to be seen
whether this is robust to additional data sets (and in the future it should be
tested against the full KiDS and BOSS likelihoods).

Before including external data, we comment on the possibility for non-thermal
hot dark matter to resolve the Hubble tension (see e.g. [55, 56, 57] for a
review). We find that, whether we include the $S_{8}$ prior or not, the value
of $H_{0}$ is barely affected by the extra $\Delta N_{\rm eff}$ (in fact, even
shifted slightly towards lower $H_{0}$ due to the well-known anti-correlation
with $m_{\rm sp}^{\rm eff}$ [31]). We therefore confirm that these models
cannot be responsible for the high-$H_{0}$ measured with some of the local
probes.

####  III.2.2 Planck+BOSS+SN1a

When the BAO/FS and SN1a data are added to the analysis, the constraints on
the thermal neutrino mass and non-thermal hot dark matter mass strengthen. We
find $m_{\rm sp}^{\rm eff}<0.73$ eV and $\Delta N_{\rm eff}<0.12$ in the
$\nu_{\rm NT}\Lambda{\rm CDM}$ model, while we get $m_{\nu}<0.04$ eV and
$\Delta N_{\rm eff}<0.27$ in the thermal case. Still, the reconstructed
$S_{8}$ value $S_{8}(\nu_{\rm NT}\Lambda{\rm CDM})=0.814_{-0.014}^{+0.017}$
and $S_{8}(\nu\Lambda{\rm CDM})=0.83\pm 0.011$ are slightly smaller than in
the Planck-only analysis. This is because the reconstructed value of
$\omega_{\rm cdm}$ is slightly smaller in the combined analysis with BAO/FS
and SN1a data, regardless of the model.

Once the prior on $S_{8}$ is added to the analysis, we again find a mild
detection of $m_{\rm eff}=0.48_{-0.36}^{+0.17}$ eV and $\Delta N_{\rm
eff}=0.0457_{-0.031}^{+0.0038}$ . However, the mean value has decreased by
$0.5\sigma$ due to the inclusion of BAO/FS and SN1a data. This reflects in a
slightly larger reconstructed $S_{8}$ value, $S_{8}(\nu_{\rm NT}\Lambda{\rm
CDM})=0.795_{-0.013}^{+0.015}$. A similar pattern is observed in the thermal
case, for which the relaxation of the constraint to $m_{\nu}<0.057$ eV is much
milder than without BAO/FS and SN1a data, while the reconstructed
$S_{8}(\nu\Lambda{\rm CDM})=0.814\pm 0.01$ is stable. Looking at
$\chi^{2}_{\rm min}$, one can see that the non-thermal case still provides a
better fit $\Delta\chi^{2}_{\rm min}(\nu_{\rm NT}\Lambda{\rm CDM})=-3.7$ than
the thermal case $\Delta\chi^{2}_{\rm min}(\nu\Lambda{\rm CDM})=-1.8$.
However, the inclusion of the $S_{8}$ prior as increased the total
$\chi^{2}_{\rm min}$ by $+4.8$ in the non-thermal case and $+6.9$ in the
thermal case. We conclude that, while Planck and $S_{8}$ data can be
accommodated in the non-thermal scenario (not in the thermal one), the BAO/FS
and SN1a data poses a serious challenge to this model.

Figure 4: Residuals of the CMB TT, EE, lensing (left panel) and matter (right
panel) power spectra with respect to $\Lambda$CDM in the bestfit
$\nu\Lambda$CDM and $\nu_{\rm NT}\Lambda$CDM models for two different datasets
(see legend). The ‘Ext’ data refers to BAO/FS+SN1a.

###  III.3 Understanding the MCMC

Figure 5: Transverse BAO (top left panel), longitudinal BAO (top right panel)
and growth factor (bottom panel) in the bestfit $\nu\Lambda$CDM and $\nu_{\rm
NT}\Lambda$CDM models for two different datasets (see legend). The ‘Ext’ data
refers to BAO/FS+SN1a. The transverse BAO has been normalized to the
$\Lambda$CDM prediction, as in Ref. [58].

In order to understand better the results of the MCMC analyses, we show in
Fig. 4 the residuals of the CMB TT, EE, lensing (top panel) and matter (bottom
panel) power spectra with respect to $\Lambda$CDM in the bestfit
$\nu\Lambda$CDM and $\nu_{\rm NT}\Lambda$CDM models obtained when considering
Planck+S8 and Planck+Ext+S8 data. We also show in Fig. 5 the corresponding
transverse BAO (top panel), longitudinal BAO (middle panel) and growth factor
(bottom panel). The first thing to notice is that, for a similar effect in the
CMB power spectra, the corresponding power suppression in the matter power
spectrum is much stronger in the $\nu_{\rm NT}\Lambda$CDM than in the
$\nu\Lambda$CDM model. This is the reason why the $\nu_{\rm NT}\Lambda$CDM can
perform much better in resolving the $S_{8}$ tension.

Looking at the BAO and $f\sigma_{8}$ prediction, one can see that the most
important difference is in the latter, which is significantly lower at all $z$
in the $\nu_{\rm NT}\Lambda$CDM because of this power suppression. This
explains the small degradation in $\chi^{2}$ in the combined analysis with
$S_{8}$. Moreover, the reconstructed dark matter density $\omega_{\rm cdm}$ in
the $\nu_{\rm NT}\Lambda$CDM is also shifted by roughly $\sim 1\sigma$
downwards (to compensate for the higher energy density due to the non-
relativistic transition of the non-thermal neutrinos), which also leads to a
small degradation in the fit to Planck data (hardly visible by eye in CMB
power spectra residuals). This small difference in the matter density is also
visible in the small$-k$ (large scales) branch of the matter power spectrum,
particularly sensitive to $\Omega_{m}$ [59]. While these differences do not
yet unambiguously rule out the $\nu_{\rm NT}\Lambda$CDM as a resolution to the
$S_{8}$ tension, they do provide an interesting avenue to probe the model with
future data, in particular through accurate measurements of the matter power
spectrum, CMB lensing power spectrum and growth factor $f\sigma_{8}$. An
interesting avenue to improve over the $\nu_{\rm NT}\Lambda$CDM results
presented here is to assume that the hot component comes from the decay of a
meta-stable cold dark matter species in the late-universe [19, 20], instead of
being present at all times. A good fit to all data can then be obtained when
the mass-ratio of the mother and daughter particle $\varepsilon\sim 0.007$ and
the CDM lifetime $\tau\sim~{}55$ Gyrs.

##  IV Implications for other non-thermal hot dark matter models

Model | Non-thermal | Thermal | Dodelson Widrow
---|---|---|---
Data set |  $m_{\rm sp}$ [eV] | $B_{\rm sp}$ |  $m_{\rm sp}$ [eV] | $\frac{T_{s}}{T_{\nu}}$ |  $m_{\rm sp}$ [eV] | $\chi$
Planck | 0.05 | 0.01 | 0 | 0.40 | 0 | 0.03
Planck+$S_{8}$ | 38.62 | 0.012 | 11.36 | 0.43 | 26.43 | 0.03
Planck+Ext | 18.98 | 0.01 | 04.59 | 0.36 | 12.85 | 0.02
Planck+Ext+$S_{8}$ | 39.81 | 0.01 | 11.75 | 0.43 | 27.49 | 0.03

Table 3: Best-fit values of the physical parameters in the non-thermal,
thermal and Dodelson-Widrow sterile neutrino models derived from our analyses.

Figure 6: Residuals of Matter power spectra (top panel) and $C_{l}$ TT, TE and
EE power spectra (bottom panel).

As discussed in the introduction and section II.2 , any distribution with the
same values of $\Delta N_{\rm eff}$ and $m_{\rm sp}^{\rm eff}$ as ours should
also relax the $\sigma_{8}$ tension. Our results can thus be used to extract
implications for the microscopic parameters of models which have momentum
distributions different from the one we have used. Here, we present such
results for two models:

a) Sterile particles at a different temperature from that of the Standard
Model neutrinos. In this model, sterile neutrinos follow a thermal Fermi-Dirac
Distribution.

$f(p)=\frac{1}{e^{p/T_{s}}+1}$ (13)

where $T_{s}$ is the temperature of sterile particles. For a thermal sterile
particle with a Fermi-Dirac distribution and a different temperature $T_{s}$,
the quantities $\Delta N_{\rm eff}$ and $\omega_{s}$ become

$\Delta N_{\rm eff}=\left(\frac{T_{s}}{T_{\nu}^{\rm
id}}\right)^{4},\qquad\omega_{s}=\frac{m_{sp}}{94.05}\left(\frac{T_{s}}{T_{\nu}}\right)^{3}\,.$
(14)

b) The Dodenson-Widrow distribution [34]

$f(p)={\chi\over{1+e^{p/T_{\nu}}}}$

where $T_{\nu}$ is the temperature of the neutrinos today, $\chi$ is a
parameter related to the phenomenological parameters as [29]:

$\Delta N_{\rm eff}=\chi,\qquad m_{\rm sp}^{\rm eff}=m_{\rm
sp}\times{\chi}\,,$ (15)

and $m_{\rm sp}$ is the individual neutrino mass in the model. We report the
best-fit value of the model parameters in Tab. 3, obtained from translating
our constraints on $\Delta N_{\rm eff}$ and $m_{\rm sp}^{\rm eff}$. We also
show in Fig. 6 the residuals of the CMB TT, TE, EE and matter power spectra
between our best-fit non-thermal HDM model and these two models. This
explicitly demonstrates our claim that, once $\Delta N_{\rm eff}$ and $m_{\rm
sp}^{\rm eff}$ are fixed, observables are indistinguishable. We note that the
residuals between the thermal neutrino model at different temperature and our
non-thermal HDM model are of the order of the sensitivity of future LSS
experiments such as EUCLID and LSST, and therefore this simple mapping might
become limited in the future. Note that, to avoid biasing constraints due to
prior effects, we refrain from translating our reconstructed posterior on
$\Delta N_{\rm eff}$ and $m_{\rm sp}^{\rm eff}$ into the model parameters.

The values we report in Tab. 3 have direct implication for thermalized hidden
sector from both particle physics [60] and cosmological perspective [61, 62].
Interestingly the main parameter for building a thermal hidden sector model is
the temperature ratio $\xi=\frac{T_{s}}{T_{vis}}$ which received a competitive
constraint (though it depends on the model) from our analysis and it may have
strong implications for light sterile neutrino [61] or other hidden sector
particle physics models [63, 64]. If the thermal hidden particle interacts
with dark matter or other particle in dark sector, the coupling and other
particle physics parameters can be constrained from our result [65].

It is tantalizing to connect the hot dark matter discussed here to the
longstanding (and debated) short base line (SBL) anomalies [66, 67] (see [68,
69] for recent reviews). Concretely, within the so called “3+1” neutrino
scenario, those can be explained by a sterile neutrino with
$m_{s}\simeq\sqrt{\Delta m_{41}^{2}}~{}1$eV and a mixing angle leading to
$\Delta N_{\rm eff}\simeq 1$. However, we find that the sterile particles
required by the $S_{8}$-tension hints to a somewhat higher mass range
$m_{s}\sim{\cal O}(10)$ eV (see tab. 3), and an almost negligible $\Delta
N_{\rm eff}$. Our constraints, whether we include the $S_{8}$-prior or not,
thus further confirm that a viable sterile neutrino solutions to the SBL
anomalies would require some additional mechanism to prohibit large $\Delta
N_{\rm eff}$ production (see e.g. [70, 71, 72, 73] for examples).
Nevertheless, it could be interesting to perform analysis including results
short base line neutrino oscillation (e.g. with an additional prior as in Ref.
[61]). This is beyond the scope of this paper and is kept for future study.

Finally, we also note that including data from the Bicep2/Kek array [74, 75],
SPT-3G [76] or ACT [77] could help further constraining the sterile neutrino
parameters thanks to higher accuracy measurement of the CMB damping tail and
lensing spectrum. We also keep that for a future study, but refer to Refs.
[78, 76] for examples (constraints typically increases by $\sim 10\%$, without
considering a prior on $S_{8}$).

##  V Conclusions

In this paper, we have explored the possibility that the ‘$S_{8}$-tension’,
the long-standing discrepancy between the determination of the amplitude of
the matter fluctuations from local [2, 3, 4, 5, 6, 7] and cosmological [1]
probes, is due to the existence of a non-thermal HDM contributing to a
fraction of the DM density in the universe and leading to a power suppression
at small-scales in the matter power spectrum. Concretely, we have considered
non-thermal HDM produced as decay products of the inflaton. Such particles
have the momentum distribution associated with decays taking place in a matter
dominated universe evolving to radiation domination, as shown in [28].
However, we have argued that any model leading to the same $\Delta N_{\rm
eff}$ and $m_{\rm sp}^{\rm eff}$ as our model (barring additional new physics
ingredients) would lead to similar effects on cosmological observables, and
therefore our constraints generically apply to any HDM models.

We have performed a comprehensive monte-carlo markov chains (MCMC) analysis
against up-to-date data from Planck, BOSS (BAO and $f\sigma_{8}$) and Pantheon
data, with and without the inclusion of a prior on the value of $S_{8}$ as
measured with the KiDS/Viking+BOSS+2dFLens data. Our findings can be
summarized as follows:

  1. •

the $\nu_{\rm NT}\Lambda$CDM model can indeed alleviate the tension between
Planck and $S_{8}$ measurements, but the success of the resolution is degraded
once BOSS and Pantheon data are included in the analysis.

  2. •

Compared to standard thermal neutrinos, the $\nu_{\rm NT}\Lambda$CDM leads to
a much stronger suppression in the matter power spectrum at late-times for
similar effect on the CMB power spectrum, and therefore to a more significant
decrease in $\sigma_{8}$.

  3. •

The impact of the $\nu_{\rm NT}\Lambda$CDM is barely visible on the BAO scale
and luminosity distance, but does significantly affect $f\sigma_{8}$
predictions. The model is therefore constrained by BOSS growth factor
measurements, and future measurements of the matter power-spectrum and
$f\sigma_{8}$ at late-times will further test this scenario.

  4. •

We further discussed the connection between our model and generic
phenomenological parameters constrained by the data that can be easily used to
translate our constraints onto other similar models. Especially, we put
constraints on other non-thermal HDM models–like the Dodelson-Widrow models or
on a thermal sterile particle with a different temperature in hidden sector.
We report competitive constraint on the hidden sector temperature and DW
scaling parameter which can have interesting particle physics implications,
for instance in the context of SBL anomalies [66, 67, 68, 69].

Lyman-$\alpha$ forest flux power spectrum data [79, 80, 81, 82, 10, 83, 84,
85] and future high accuracy measurement of the matter power spectrum at small
scales by upcoming surveys such as Euclid [86], LSST [87], and DESI [88] can
further test these models as a resolution to the $S_{8}$-tension.

## acknowledgements

We thank Guillermo Franco Abellàn and Riccardo Murgia for many interesting
discussions. AM is supported in part by the SERB, DST, Government of India by
the grant MTR/2019/000267. SD acknowledges SERB DST Government of India grant
CRG/2019/006147 for supporting the project. VP is partly supported by the
CNRS-IN2P3 grant Dark21 and by the European Union’s Horizon 2020 research and
innovation program under the Marie Skodowska-Curie grant agreement No
860881-HIDDeN. The authors acknowledge the use of computational resources from
the Dark Energy computing Center funded by the OCEVU Labex (ANR-11-LABX-0060)
and the Excellence Initiative of Aix-Marseille University (A*MIDEX) of the
“Investissements d’Avenir” programme as well as IIA Nova cluster where initial
analysis was carried out.

##  Appendix A $\chi^{2}_{\rm min}$ per experiment

We report $\chi^{2}_{\rm min}$ per experiment in each of the analysis
performed.

Experiment |  $\Lambda$CDM
---|---
Planck high$-\ell$ TT,TE,EE | 2346.7 | 2350.8 | 2346 | 2349.1
Planck low$-\ell$ EE | 396 | 396.1 | 396.8 | 396.2
Planck low$-\ell$ TT | 23.2 | 22.5 | 23.4 | 22.6
Planck lensing | 8.8 | 9.6 | 9.2 | 9.1
Pantheon | $-$ | $-$ | 1026.9 | 1026.7
BAO/FS BOSS DR12 | $-$ | $-$ | 6.9 | 6.5
BAO BOSS low$-z$ | $-$ | $-$ | 1.2 | 2.3
KiDS/BOSS/2dFGS | $-$ | 4.6 | $-$ | 5.9
total | 2774.8 | 2783.4 | 3810.4 | 3818.2

Table 4: Best-fit $\chi^{2}$ per experiment (and total) in the $\Lambda$CDM
model.

Experiment |  $\nu\Lambda$CDM
---|---
Planck high$-\ell$ TT,TE,EE | 2345.98 | 2348.2 | 2346.9 | 2348.6
Planck low$-\ell$ EE | 396.54 | 396.8 | 396.5 | 395.7
Planck low$-\ell$ TT | 23.3 | 22.2 | 22.8 | 22.4
Planck lensing | 9.03 | 8.9 | 8.8 | 9.3
Pantheon | $-$ | $-$ | 1026.8 | 1026.7
BAO/FS BOSS DR12 | $-$ | $-$ | 6.1 | 5.9
BAO BOSS low$-z$ | $-$ | $-$ | 1.7 | 1.7
KiDS/BOSS/2dFGS | $-$ | 5.8 | $-$ | 6.1
total | 2774.9 | 2782.0 | 3809.5 | 3816.4

Table 5: Best-fit $\chi^{2}$ per experiment (and total) in the model with
massive thermal neutrinos and additional relativistic degrees of freedom.

Experiment |  $\nu_{\rm NT}\Lambda$CDM
---|---
Planck high$-\ell$ TT,TE,EE | 2346.7 | 2348.7 | 2 2346.4 | 2349.1
Planck low$-\ell$ EE | 396.3 | 395.9 | 396.8 | 396.9
Planck low$-\ell$ TT | 23.1 | 23.3 | 23.4 | 23.1
Planck lensing | 8.8 | 9.2 | 8.8 | 9.1
Pantheon | $-$ | $-$ | 1026.8 | 1026.7
BAO/FS BOSS DR12 | $-$ | $-$ | 6.1 | 6.8
BAO BOSS low$-z$ | $-$ | $-$ | 1.4 | 1.7
KiDS/BOSS/2dFGS | $-$ | 1.6 | $-$ | 1.2
total | 2775.0 | 2778.6 | 3809.7 | 3814.5

Table 6: Best-fit $\chi^{2}$ per experiment (and total) in the non-thermal
sterile neutrino model.

## References

  * Aghanim  _et al._ [2018a] N. Aghanim _et al._ (Planck), (2018a), 1807.06209 .
  * Heymans  _et al._ [2013] C. Heymans _et al._ , CFHTLenS tomographic weak lensing cosmological parameter constraints: Mitigating the impact of intrinsic galaxy alignments, Mon. Not. Roy. Astron. Soc. 432, 2433 (2013), arXiv:1303.1808 [astro-ph.CO] .
  * MacCrann  _et al._ [2015] N. MacCrann, J. Zuntz, S. Bridle, B. Jain, and M. R. Becker, Cosmic Discordance: Are Planck CMB and CFHTLenS weak lensing measurements out of tune?, Mon. Not. Roy. Astron. Soc. 451, 2877 (2015), arXiv:1408.4742 [astro-ph.CO] .
  * Hildebrandt  _et al._ [2020] H. Hildebrandt _et al._ , KiDS+VIKING-450: Cosmic shear tomography with optical and infrared data, Astron. Astrophys. 633, A69 (2020), arXiv:1812.06076 [astro-ph.CO] .
  * Joudaki  _et al._ [2020] S. Joudaki _et al._ , KiDS+VIKING-450 and DES-Y1 combined: Cosmology with cosmic shear, Astron. Astrophys. 638, L1 (2020), arXiv:1906.09262 [astro-ph.CO] .
  * Abbott  _et al._ [2018] T. Abbott _et al._ (DES), Dark Energy Survey year 1 results: Cosmological constraints from galaxy clustering and weak lensing, Phys. Rev. D 98, 043526 (2018), arXiv:1708.01530 [astro-ph.CO] .
  * Heymans  _et al._ [2020] C. Heymans _et al._ , KiDS-1000 Cosmology: Multi-probe weak gravitational lensing and spectroscopic galaxy clustering constraints (2020) arXiv:2007.15632 [astro-ph.CO] .
  * Kumar and Nunes [2016] S. Kumar and R. C. Nunes, Phys. Rev. D94, 123511 (2016), 1608.02454 .
  * Murgia  _et al._ [2016] R. Murgia, S. Gariazzo, and N. Fornengo, JCAP 04, 014.
  * Archidiacono  _et al._ [2019] M. Archidiacono, D. C. Hooper, R. Murgia, S. Bohr, J. Lesgourgues, and M. Viel, Constraining Dark Matter-Dark Radiation interactions with CMB, BAO, and Lyman-$\alpha$, JCAP 10, 055, arXiv:1907.01496 [astro-ph.CO] .
  * Di Valentino  _et al._ [2020] E. Di Valentino _et al._ , (2020), 2008.11285 .
  * Becker  _et al._ [2021] N. Becker, D. C. Hooper, F. Kahlhoefer, J. Lesgourgues, and N. Schöneberg, Cosmological constraints on multi-interacting dark matter, JCAP 02, 019, arXiv:2010.04074 [astro-ph.CO] .
  * Enqvist  _et al._ [2015] K. Enqvist, S. Nadathur, T. Sekiguchi, and T. Takahashi, JCAP 09, 067, 1505.05511 .
  * Poulin  _et al._ [2016] V. Poulin, P. D. Serpico, and J. Lesgourgues, JCAP 1608 (08), 036, 1606.02073 .
  * Vattis  _et al._ [2019] K. Vattis, S. M. Koushiappas, and A. Loeb, Phys. Rev. D99, 121302 (2019), 1903.06220 .
  * Haridasu and Viel [2020] B. S. Haridasu and M. Viel, (2020), 2004.07709 .
  * Clark  _et al._ [2020] S. J. Clark, K. Vattis, and S. M. Koushiappas, (2020), 2006.03678 .
  * Pandey  _et al._ [2020] K. L. Pandey, T. Karwal, and S. Das, Alleviating the $H_{0}$ and $\sigma_{8}$ anomalies with a decaying dark matter model, JCAP 07, 026, arXiv:1902.10636 [astro-ph.CO] .
  * Abellan  _et al._ [2020] G. F. Abellan, R. Murgia, V. Poulin, and J. Lavalle, (2020), 2008.09615 .
  * Abellán  _et al._ [2021] G. F. Abellán, R. Murgia, and V. Poulin, Linear cosmological constraints on 2-body decaying dark matter scenarios and robustness of the resolution to the $S_{8}$ tension, (2021), arXiv:2102.12498 [astro-ph.CO] .
  * Poulin  _et al._ [2018] V. Poulin, K. K. Boddy, S. Bird, and M. Kamionkowski, Implications of an extended dark energy cosmology with massive neutrinos for cosmological tensions, Phys. Rev. D97, 123504 (2018), arXiv:1803.02474 [astro-ph.CO] .
  * Kreisch  _et al._ [2019] C. D. Kreisch, F.-Y. Cyr-Racine, and O. Doré, The Neutrino Puzzle: Anomalies, Interactions, and Cosmological Tensions, (2019), arXiv:1902.00534 [astro-ph.CO] .
  * Boyarsky  _et al._ [2021] A. Boyarsky, M. Ovchynnikov, N. Sabti, and V. Syvolap, When FIMPs Decay into Neutrinos: The $N_{\mathrm{eff}}$ Story, (2021), arXiv:2103.09831 [hep-ph] .
  * Hou  _et al._ [2013] Z. Hou, R. Keisler, L. Knox, M. Millea, and C. Reichardt, How Massless Neutrinos Affect the Cosmic Microwave Background Damping Tail, Phys. Rev. D 87, 083008 (2013), arXiv:1104.2333 [astro-ph.CO] .
  * Scherrer and Turner [1988] R. J. Scherrer and M. S. Turner, Primordial Nucleosynthesis with Decaying Particles. 1. Entropy Producing Decays. 2. Inert Decays, Astrophys. J. 331, 19 (1988).
  * Hasenkamp and Kersten [2013] J. Hasenkamp and J. Kersten, Dark radiation from particle decay: cosmological constraints and opportunities, JCAP 08, 024, arXiv:1212.4160 [hep-ph] .
  * Conlon and Marsh [2013] J. P. Conlon and M. C. D. Marsh, The Cosmophenomenology of Axionic Dark Radiation, JHEP 10, 214, arXiv:1304.1804 [hep-ph] .
  * Bhattacharya  _et al._ [2020] S. Bhattacharya, S. Das, K. Dutta, M. R. Gangopadhyay, R. Mahanta, and A. Maharana, Non-thermal Hot Dark Matter from Inflaton/Moduli Decay: The Momentum Distribution and Relaxing the Cosmological Mass Bound, (2020), arXiv:2009.05987 [astro-ph.CO] .
  * Acero and Lesgourgues [2009] M. A. Acero and J. Lesgourgues, Cosmological constraints on a light non-thermal sterile neutrino, Phys. Rev. D 79, 045026 (2009), arXiv:0812.2249 [astro-ph] .
  * Lesgourgues and Pastor [2006] J. Lesgourgues and S. Pastor, Massive neutrinos and cosmology, Phys. Rept. 429, 307 (2006), astro-ph/0603494 .
  * Lesgourgues [2011] J. Lesgourgues, The Cosmic Linear Anisotropy Solving System (CLASS) I: Overview, (2011), arXiv:1104.2932 [astro-ph.IM] .
  * Roy Choudhury and Hannestad [2020] S. Roy Choudhury and S. Hannestad, Updated results on neutrino mass and mass hierarchy from cosmology with Planck 2018 likelihoods, JCAP 07, 037, arXiv:1907.12598 [astro-ph.CO] .
  * Roy Choudhury and Choubey [2019a] S. Roy Choudhury and S. Choubey, Constraining light sterile neutrino mass with the BICEP2/Keck Array 2014 B-mode polarization data, Eur. Phys. J. C 79, 557 (2019a), arXiv:1807.10294 [astro-ph.CO] .
  * Dodelson and Widrow [1994] S. Dodelson and L. M. Widrow, Sterile neutrinos as dark matter, Physical Review Letters 72, 17 (1994).
  * Cuoco  _et al._ [2005] A. Cuoco, J. Lesgourgues, G. Mangano, and S. Pastor, Do observations prove that cosmological neutrinos are thermally distributed?, Phys. Rev. D 71, 123501 (2005), arXiv:astro-ph/0502465 .
  * Feng  _et al._ [2008] J. L. Feng, H. Tu, and H.-B. Yu, Thermal Relics in Hidden Sectors, JCAP 10, 043, arXiv:0808.2318 [hep-ph] .
  * Das and Sigurdson [2012] S. Das and K. Sigurdson, Cosmological Limits on Hidden Sector Dark Matter, Phys. Rev. D 85, 063510 (2012), arXiv:1012.4458 [astro-ph.CO] .
  * Berlin  _et al._ [2019] A. Berlin, N. Blinov, and S. W. Li, Dark Sector Equilibration During Nucleosynthesis, Phys. Rev. D 100, 015038 (2019), arXiv:1904.04256 [hep-ph] .
  * Gelmini  _et al._ [2004] G. Gelmini, S. Palomares-Ruiz, and S. Pascoli, Low reheating temperature and the visible sterile neutrino, Phys. Rev. Lett. 93, 081302 (2004), arXiv:astro-ph/0403323 .
  * Gelmini  _et al._ [2008] G. Gelmini, E. Osoba, S. Palomares-Ruiz, and S. Pascoli, MeV sterile neutrinos in low reheating temperature cosmological scenarios, JCAP 10, 029, arXiv:0803.2735 [astro-ph] .
  * Lesgourgues [2011] J. Lesgourgues, The Cosmic Linear Anisotropy Solving System (CLASS) I: Overview, arXiv e-prints , arXiv:1104.2932 (2011), arXiv:1104.2932 [astro-ph.IM] .
  * Lesgourgues and Tram [2011] J. Lesgourgues and T. Tram, The Cosmic Linear Anisotropy Solving System (CLASS) IV: efficient implementation of non-cold relics, JCAP 2011 (9), 032, arXiv:1104.2935 [astro-ph.CO] .
  * Lesgourgues and Tram [2011] J. Lesgourgues and T. Tram, The Cosmic Linear Anisotropy Solving System (CLASS) IV: efficient implementation of non-cold relics, JCAP 09, 032, arXiv:1104.2935 [astro-ph.CO] .
  * Aghanim  _et al._ [2018b] N. Aghanim _et al._ (Planck), Planck 2018 results. VI. Cosmological parameters, (2018b), arXiv:1807.06209 [astro-ph.CO] .
  * Beutler  _et al._ [2011] F. Beutler, C. Blake, M. Colless, D. H. Jones, L. Staveley-Smith, L. Campbell, Q. Parker, W. Saunders, and F. Watson, The 6dF Galaxy Survey: Baryon Acoustic Oscillations and the Local Hubble Constant, Mon. Not. Roy. Astron. Soc. 416, 3017 (2011), arXiv:1106.3366 [astro-ph.CO] .
  * Ross  _et al._ [2015] A. J. Ross, L. Samushia, C. Howlett, W. J. Percival, A. Burden, and M. Manera, The clustering of the SDSS DR7 main Galaxy sample – I. A 4 per cent distance measure at $z=0.15$, Mon. Not. Roy. Astron. Soc. 449, 835 (2015), arXiv:1409.3242 [astro-ph.CO] .
  * Alam  _et al._ [2017] S. Alam _et al._ (BOSS), The clustering of galaxies in the completed SDSS-III Baryon Oscillation Spectroscopic Survey: cosmological analysis of the DR12 galaxy sample, Mon. Not. Roy. Astron. Soc. 470, 2617 (2017), arXiv:1607.03155 [astro-ph.CO] .
  * de Sainte Agathe  _et al._ [2019] V. de Sainte Agathe _et al._ , Baryon acoustic oscillations at z = 2.34 from the correlations of Ly$\alpha$ absorption in eBOSS DR14, Astron. Astrophys. 629, A85 (2019), arXiv:1904.03400 [astro-ph.CO] .
  * Blomqvist  _et al._ [2019] M. Blomqvist _et al._ , Baryon acoustic oscillations from the cross-correlation of Ly$\alpha$ absorption and quasars in eBOSS DR14, Astron. Astrophys. 629, A86 (2019), arXiv:1904.03430 [astro-ph.CO] .
  * Scolnic  _et al._ [2018] D. M. Scolnic _et al._ , The Complete Light-curve Sample of Spectroscopically Confirmed SNe Ia from Pan-STARRS1 and Cosmological Constraints from the Combined Pantheon Sample, Astrophys. J. 859, 101 (2018), arXiv:1710.00845 [astro-ph.CO] .
  * Brinckmann and Lesgourgues [2018] T. Brinckmann and J. Lesgourgues, MontePython 3: boosted MCMC sampler and other features, (2018), arXiv:1804.07261 [astro-ph.CO] .
  * James and Roos [1975] F. James and M. Roos, Minuit – a system for function minimization and analysis of the parameter errors and correlations, Computer Physics Communications 10, 343 (1975).
  * Lewis  _et al._ [2000] A. Lewis, A. Challinor, and A. Lasenby, Efficient computation of CMB anisotropies in closed FRW models, Astrophys. J. 538, 473 (2000), arXiv:astro-ph/9911177 [astro-ph] .
  * Gelman and Rubin [1992] A. Gelman and D. B. Rubin, Inference from Iterative Simulation Using Multiple Sequences, Statist. Sci. 7, 457 (1992).
  * Riess  _et al._ [2020] A. G. Riess, S. Casertano, W. Yuan, J. B. Bowers, L. Macri, J. C. Zinn, and D. Scolnic, (2020), 2012.08534 .
  * Freedman  _et al._ [2020] W. L. Freedman, B. F. Madore, T. Hoyt, I. S. Jang, R. Beaton, M. G. Lee, A. Monson, J. Neeley, and J. Rich, Calibration of the Tip of the Red Giant Branch (TRGB), (2020), arXiv:2002.01550 [astro-ph.GA] .
  * Di Valentino  _et al._ [2021] E. Di Valentino, O. Mena, S. Pan, L. Visinelli, W. Yang, A. Melchiorri, D. F. Mota, A. G. Riess, and J. Silk, In the Realm of the Hubble tension $-$ a Review of Solutions, (2021), arXiv:2103.01183 [astro-ph.CO] .
  * Aghanim  _et al._ [2020] N. Aghanim _et al._ (Planck), Planck 2018 results. I. Overview and the cosmological legacy of Planck, Astron. Astrophys. 641, A1 (2020), arXiv:1807.06205 [astro-ph.CO] .
  * Lesgourgues  _et al._ [2013] J. Lesgourgues, G. Mangano, G. Miele, and S. Pastor,  _Neutrino Cosmology_ (Cambridge University Press, 2013).
  * Reece and Roxlo [2016] M. Reece and T. Roxlo, Nonthermal production of dark radiation and dark matter, JHEP 09, 096, arXiv:1511.06768 [hep-ph] .
  * Gariazzo  _et al._ [2013] S. Gariazzo, C. Giunti, and M. Laveder, Light Sterile Neutrinos in Cosmology and Short-Baseline Oscillation Experiments, JHEP 11, 211, arXiv:1309.3192 [hep-ph] .
  * Cyr-Racine  _et al._ [2014] F.-Y. Cyr-Racine, R. de Putter, A. Raccanelli, and K. Sigurdson, Constraints on Large-Scale Dark Acoustic Oscillations from Cosmology, Phys. Rev. D 89, 063517 (2014), arXiv:1310.3278 [astro-ph.CO] .
  * Foot [2014] R. Foot, Mirror dark matter: Cosmology, galaxy structure and direct detection, Int. J. Mod. Phys. A 29, 1430013 (2014), arXiv:1401.3965 [astro-ph.CO] .
  * Franχ_bca  _et al._ [2013] U. Franχ_bca, R. A. Lineros, J. Palacio, and S. Pastor, Probing interactions within the dark matter sector via extra radiation contributions, Phys. Rev. D 87, 123521 (2013), arXiv:1303.1776 [astro-ph.CO] .
  * Brust  _et al._ [2013] C. Brust, D. E. Kaplan, and M. T. Walters, New Light Species and the CMB, JHEP 12, 058, arXiv:1303.5379 [hep-ph] .
  * Aguilar-Arevalo  _et al._ [2001] A. Aguilar-Arevalo _et al._ (LSND), Evidence for neutrino oscillations from the observation of $\bar{\nu}_{e}$ appearance in a $\bar{\nu}_{\mu}$ beam, Phys. Rev. D 64, 112007 (2001), arXiv:hep-ex/0104049 .
  * Aguilar-Arevalo  _et al._ [2013] A. A. Aguilar-Arevalo _et al._ (MiniBooNE), Improved Search for $\bar{\nu}_{\mu}\rightarrow\bar{\nu}_{e}$ Oscillations in the MiniBooNE Experiment, Phys. Rev. Lett. 110, 161801 (2013), arXiv:1303.2588 [hep-ex] .
  * [68] M. Maltoni.
  * Dentler  _et al._ [2018] M. Dentler, A. Hernández-Cabezudo, J. Kopp, P. A. N. Machado, M. Maltoni, I. Martinez-Soler, and T. Schwetz, Updated Global Analysis of Neutrino Oscillations in the Presence of eV-Scale Sterile Neutrinos, JHEP 08, 010, arXiv:1803.10661 [hep-ph] .
  * Hamann  _et al._ [2011] J. Hamann, S. Hannestad, G. G. Raffelt, and Y. Y. Y. Wong, Sterile neutrinos with eV masses in cosmology: How disfavoured exactly?, JCAP 09, 034, arXiv:1108.4136 [astro-ph.CO] .
  * Archidiacono  _et al._ [2016] M. Archidiacono, S. Gariazzo, C. Giunti, S. Hannestad, R. Hansen, M. Laveder, and T. Tram, Pseudoscalar—sterile neutrino interactions: reconciling the cosmos with neutrino oscillations, JCAP 08, 067, arXiv:1606.07673 [astro-ph.CO] .
  * Chu  _et al._ [2015] X. Chu, B. Dasgupta, and J. Kopp, Sterile neutrinos with secret interactions—lasting friendship with cosmology, JCAP 10, 011, arXiv:1505.02795 [hep-ph] .
  * de Salas  _et al._ [2015] P. F. de Salas, M. Lattanzi, G. Mangano, G. Miele, S. Pastor, and O. Pisanti, Bounds on very low reheating scenarios after Planck, Phys. Rev. D 92, 123534 (2015), arXiv:1511.00672 [astro-ph.CO] .
  * Ade  _et al._ [2016] P. A. R. Ade _et al._ (BICEP2, Keck Array), Improved Constraints on Cosmology and Foregrounds from BICEP2 and Keck Array Cosmic Microwave Background Data with Inclusion of 95 GHz Band, Phys. Rev. Lett. 116, 031302 (2016), arXiv:1510.09217 [astro-ph.CO] .
  * Ade  _et al._ [2018] P. A. R. Ade _et al._ (BICEP2, Keck Array), BICEP2 / Keck Array x: Constraints on Primordial Gravitational Waves using Planck, WMAP, and New BICEP2/Keck Observations through the 2015 Season, Phys. Rev. Lett. 121, 221301 (2018), arXiv:1810.05216 [astro-ph.CO] .
  * Balkenhol  _et al._ [2021] L. Balkenhol _et al._ (SPT), Constraints on $\Lambda$CDM Extensions from the SPT-3G 2018 $EE$ and $TE$ Power Spectra, (2021), arXiv:2103.13618 [astro-ph.CO] .
  * Aiola  _et al._ [2020] S. Aiola _et al._ (ACT), The Atacama Cosmology Telescope: DR4 Maps and Cosmological Parameters, (2020), arXiv:2007.07288 [astro-ph.CO] .
  * Roy Choudhury and Choubey [2019b] S. Roy Choudhury and S. Choubey, Constraining light sterile neutrino mass with the BICEP2/Keck Array 2014 B-mode polarization data, Eur. Phys. J. C 79, 557 (2019b), arXiv:1807.10294 [astro-ph.CO] .
  * Wang  _et al._ [2013] M.-Y. Wang, R. A. C. Croft, A. H. G. Peter, A. R. Zentner, and C. W. Purcell, Phys. Rev. D 88, 123515 (2013), 1309.7354 .
  * Murgia  _et al._ [2017] R. Murgia, A. Merle, M. Viel, M. Totzauer, and A. Schneider, ”Non-cold” dark matter at small scales: a general approach, JCAP 11, 046, arXiv:1704.07838 [astro-ph.CO] .
  * Baur  _et al._ [2017] J. Baur, N. Palanque-Delabrouille, C. Yeche, A. Boyarsky, O. Ruchayskiy, E. Armengaud, and J. Lesgourgues, Constraints from Ly-$\alpha$ forests on non-thermal dark matter including resonantly-produced sterile neutrinos, JCAP 12, 013, arXiv:1706.03118 [astro-ph.CO] .
  * Murgia  _et al._ [2018] R. Murgia, V. Irsic, and M. Viel, Novel constraints on noncold, nonthermal dark matter from Lyman- $\alpha$ forest data, Phys. Rev. D 98, 083540 (2018), arXiv:1806.08371 [astro-ph.CO] .
  * Miller  _et al._ [2019] C. Miller, A. L. Erickcek, and R. Murgia, Constraining nonthermal dark matter’s impact on the matter power spectrum, Phys. Rev. D 100, 123520 (2019), arXiv:1908.10369 [astro-ph.CO] .
  * Palanque-Delabrouille  _et al._ [2020] N. Palanque-Delabrouille, C. Yèche, N. Schöneberg, J. Lesgourgues, M. Walther, S. Chabanier, and E. Armengaud, Hints, neutrino bounds and WDM constraints from SDSS DR14 Lyman-$\alpha$ and Planck full-survey data, JCAP 04, 038, arXiv:1911.09073 [astro-ph.CO] .
  * Enzi  _et al._ [2020] W. Enzi _et al._ , (2020), 2010.13802 .
  * Amendola  _et al._ [2018] L. Amendola _et al._ , Cosmology and fundamental physics with the Euclid satellite, Living Rev. Rel. 21, 2 (2018), arXiv:1606.00180 [astro-ph.CO] .
  * Alonso  _et al._ [2018] D. Alonso _et al._ (LSST Dark Energy Science), The LSST Dark Energy Science Collaboration (DESC) Science Requirements Document, (2018), arXiv:1809.01669 [astro-ph.CO] .
  * Aghamousa  _et al._ [2016] A. Aghamousa _et al._ (DESI), The DESI Experiment Part I: Science,Targeting, and Survey Design, (2016), arXiv:1611.00036 [astro-ph.IM] .