$P = 4kT$ (where $k$ = Boltzmann’s constant, and $T$ = temperature of the instrument ($K$))
And the mean voltage is thus, of course, $V^2/R = 4kT$
And the voltage is distributed as a Gaussian around that mean, with RMS given by $\sqrt{4kRT\Delta f}$.
Emphasis mine:
I want to plot the expected voltage vs time trace
The expectation (that's a term from probability theory and statistics) of a white process is zero. So, your expected voltage trace is a constant 0.
You might mean that you want to plot an exemplary realization of a random noise process - in that case:
Notice that in either case, the drawing you made has probability 0 of actually being the realization of the white noise you want to be observing; you can't deterministically draw a random thing.
Before getting to the answer, several corrections are necessary. First, the units of the Boltzmann constant are J/K, so kT is energy, not power. Thus P is not equal to 4kT. Second, the expected value of the thermal noise voltage, $\mu$, is zero, as required by the second law of thermodynamics. So $V^2/R$ is not equal to 4kT.
When I do computer simulations involving Gaussian white noise, of which thermal noise is a type, I assume that I am sampling a pseudo-random Gaussian noise generator, i.e., the software I am using provides a way to get successive negligibly-correlated samples. I assume a constant sampling rate of $f_s$ samples per second. The constant point spacing, $\Delta t$, is $1/f_s$ seconds. Then the Nyquist frequency, $f_{nyq}$, is $f_s /2$ Hz, which is equivalent to $1/2\Delta t$ Hz.
The Gaussian white noise has zero population mean, $\mu$, and $\sigma^2$ population variance, so the bilateral noise power spectral density, $\eta$, is $\sigma^2 /f_s$. This simply means that the total noise power, $\sigma^2$, is uniformly distributed, because the noise is white, from minus $f_{nyq}$ to plus $f_{nyq}$. This frequency interval is simply $f_s$. Hence
$$ \eta = \sigma^2 /f_s = \sigma^2 /2f_{nyq} = \sigma^2 Δt \tag 1 $$
The unilateral noise power spectral density, $2 \eta$, is $2 \sigma^2 Δt$. Therefore
$$ 2 \eta = 2\sigma^2 /f_s = \sigma^2 /f_{nyq} = 2\sigma^2 Δt \tag 2 $$
For thermal noise below roughly 10 THz (and see here for a bit more on that), the unilateral noise power spectral density is well approximated as $4RkT$. Hence
$$ 4RkT = 2 \eta = 2\sigma^2 Δt = 2\sigma^2 /f_s = \sigma^2 /f_{nyq} \tag 3 $$
Thus
$$ \sigma^2 = 2RkT / Δt = 2RkTf_s = 4RkTf_{nyq} \tag 4 $$
So select R, T and either $f_s$ or $\Delta t$. Then equation (4) gives the population variance for the pseudo-random Gaussian noise generator.
Example: Assume $R = 10 k \Omega$, $T = 300 K$, $k = 1.380649 \times 10^{-23} J/K$ and $\Delta t = 1 \mu s$, so $f_s = 1 MHz$. Then $\sigma^2 = 8.283894 \times 10^{-11}$ $V^2$ and $\sigma = 9.101590 \times 10^{-6} \space V$. In the software I use, the syntax for generating a sample from the pseudo-random Gaussian noise generator is simply $Gaussian(\mu ,\sigma)$. Since $\mu = 0$, and with $\sigma$ as above, this becomes $Gaussian(0,9.101590 \times 10^{-6})$. A specific temporal trace, from 0 to 1.023 ms in $1\mu s$ increments, is
From equation (3), $4RkT = 1.6567788 \times 10^{-16}$ $V^2/Hz$. This is the theoretical unilateral noise power spectral density (PSD). The plot below shows the average of $10^4$ unilateral, mean-subtracted 1024 point PSDs. It is seen that the PSD is white and in excellent agreement with the theoretical value.
