3 added 26 characters in body

The answer is almost certainly that the limiting density does not exist. Without going into the details of the proof allow me to give a heuristic argument which is based on how the OP likely generated his graph of the relative frequency of happy numbers.

Let $Y_n$ be the r.v. uniformly distributed amongst integers in the interval $[0,10^n -1]$ (that is $Y_n$ picks a random $n$-digit integer). If $X_i$ denotes the r.v. for the digit of $10^i$ in $Y_n$, then $s(Y_n) = \sum\limits_{i=0}^n sum\limits_{i=0}^{n-1} s(X_i)$.

I'm guessing the way you generated your graph was you first computed the distribution of $s(Y_n)$ (this can be done recursively) then computed $\mathbb{P}\big(s(Y_n) \text{ is happy}\big)$. This would give the relative density of happy numbers amongst all $n$ digit integers.

Studying the distribution of $s(Y_n)$ can tell us a lot. Its equivalent to rolling $n$ times a 10-sided die with faces $0,1,4,\dots, 81$ and finding the sum. Its distribution is Gaussian as $n$ gets large by the central limit theorem. More importantly most of the distribution is concentrated near the mean, which is $28.5n$. This implies that the density happy numbers amongst all $n$-digit integers depends almost entirely on the distribution of happy numbers near $28.5n$.

For example, there is a peak in your graph at around $n = 400$ of about $.185$ density. Calculating the density of happy numbers within one standard deviation from the mean of $s(Y_{400})$ we get a density of .1911 (the interval I looked at was $[10916,11884]$). If you assume $s(Y_{400})$ is "exactly" normally distributed and estimated the density in this manner you would get a much better approximation.

This means picking $n$ s.t. the mean of $s(Y_n)$ lands in the interval $[10^{400},10^{401}-1]$ then the density of happy numbers amongst $n$-digit integers should be around $.185$. Likely some choices of $n$ will give densities strictly larger than $.185$ and some strictly smaller. This has led me to suspect that by iterating the this process, the upper density of happy numbers may be $1$, and lower density $0$.

The article Joe Silverman mentioned is my own. In it I attempt to give the above heuristic a rigorous foundation. It is still a rough draft and has only been reviewed by one of my fellow graduate students, so I won't to say it is definitely correct, although I am very confident it is. I have been working on it for the past few weeks, seeing your question on MO I decided to go ahead and upload a rough draft. In it I use an averaging argument to say that if you find experimentally a large interval of $n$-digit integers ($n$ sufficiently large) which contain happy numbers with density $d$, then the upper density of happy numbers is at least $d(1 - o(1))$. That is where the upper density $\geq .18$ and lower density $\leq .12$ comes from.

2 added 5 characters in body

The answer is almost certainly that the limiting density does not exist. Without going into the details of the proof allow me to give a heuristic argument which is based on how the OP likely generated his graph of the relative frequency of happy numbers.

Let $Y_n$ be the r.v. uniformly distributed amongst integers in the interval $[0,10^n -1]$ (that is $Y_n$ picks a random $n$-digit integer). If $X_i$ denotes the r.v. for the digit of $10^i$ in $Y_n$, then $s(Y_n) = \sum\limits_{i=0}^n s(X_i)$.

I'm guessing the way you generated your graph was you first computed the distribution of $s(Y_n)$ (this can be done recursively) then computed $\mathbb{P}\big(s(Y_n) \text{ is happy}\big)$. This would give the relative density of happy numbers amongst all $n$ digit integers.

Studying the distribution of $s(Y_n)$ can tell us a lot. Its equivalent to rolling $n$ times a 10-sided die with faces $0,1,4,\dots, 81$. Its distribution is Gaussian as $n$ gets large by the central limit theorem. More importantly most of the distribution is concentrated near the mean, which is $28.5n$. This implies that the density happy numbers amongst all $n$-digit integers depends almost entirely on the distribution of happy numbers near $28.5n$.

For example, there is a peak in your graph at around $n = 400$ of about $.185$ density. Calculating the density of happy numbers within one standard deviation from the mean of $s(Y_{400})$ we get a density of .1911 (the interval I looked at was $[10916,11884]$). If you assume $s(Y_{400})$ is "exactly" normally distributed and estimated the density in this manner you would get a much better approximation.

This means picking $n$ s.t. the mean of $s(n)$ s(Y_n)$lands in the interval$[10^{400},10^{401}-1]$then the density of happy numbers amongst$n$-digit integers should be around$.185$. Likely some choices of$n$will give densities strictly larger than$.185$and some strictly smaller. This has led me to suspect that by iterating the process the upper density of happy numbers may be$1$, and lower density$0$. The article Joe Silverman mentioned is my own. In it I attempt to give the above heuristic a rigorous foundation. It is still a rough draft and has only been reviewed by one of my fellow graduate students, so I won't to say it is definitely correct, although I am very confident it is. I have been working on it for the past few weeks, seeing your question on MO I decided to go ahead and upload a rough draft. In it I use an averaging argument to say that if you find experimentally a large interval of$n$-digit integers ($n$sufficiently large) which contain happy numbers with density$d$, then the upper density of happy numbers is at least$d(1 - o(1))$. That is where the upper density$\geq .18$and lower density$ \leq .12$comes from. 1 The answer is almost certainly that the limiting density does not exist. Without going into the details of the proof allow me to give a heuristic argument which is based on how the OP likely generated his graph of the relative frequency of happy numbers. Let$Y_n$be the r.v. uniformly distributed amongst integers in the interval$[0,10^n -1]$(that is$Y_n$picks a random$n$-digit integer). If$X_i$denotes the r.v. for the digit of$10^i$in$Y_n$, then$s(Y_n) = \sum\limits_{i=0}^n s(X_i)$. I'm guessing the way you generated your graph was you first computed the distribution of$s(Y_n)$(this can be done recursively) then computed$\mathbb{P}\big(s(Y_n) \text{ is happy}\big)$. This would give the relative density of happy numbers amongst all$n$digit integers. Studying the distribution of$s(Y_n)$can tell us a lot. Its equivalent to rolling$n$times a 10-sided die with faces$0,1,4,\dots, 81$. Its distribution is Gaussian as$n$gets large by the central limit theorem. More importantly most of the distribution is concentrated near the mean, which is$28.5n$. This implies that the density happy numbers amongst all$n$-digit integers depends almost entirely on the distribution of happy numbers near$28.5n$. For example, there is a peak in your graph at around$n = 400$of about$.185$density. Calculating the density of happy numbers within one standard deviation from the mean of$s(Y_{400})$we get a density of .1911 (the interval I looked at was$[10916,11884]$). If you assume$s(Y_{400})$is "exactly" normally distributed and estimated the density in this manner you would get a much better approximation. This means picking$n$s.t. the mean of$s(n)$lands in the interval$[10^{400},10^{401}-1]$then the density of happy numbers amongst$n$-digit integers should be around$.185$. Likely some choices of$n$will give densities strictly larger than$.185$and some strictly smaller. This has led me to suspect that by iterating the process the upper density happy numbers may be$1$, and lower density$0$. The article Joe Silverman mentioned is my own. In it I attempt to give the above heuristic a rigorous foundation. It is still a rough draft and has only been reviewed by one of my fellow graduate students, so I won't to say it is definitely correct, although I am very confident it is. I have been working on it for the past few weeks, seeing your question on MO I decided to go ahead and upload a rough draft. In it I use an averaging argument to say that if you find experimentally a large interval of$n$-digit integers ($n$sufficiently large) which contain happy numbers with density$d$, then the upper density of happy numbers is at least$d(1 - o(1))$. That is where the upper density$\geq .18$and lower density$ \leq .12\$ comes from.