Autocorrelation of LTI system output
Actually, we have two kinds of autocorrelation functions. One is defined for stochastic signals and the other for deterministic ones. 

If $x(t)$ is a stochastic signal, then its autocorrelation function will be $$R_x(t+\tau,t)=\Bbb E\{x(t+\tau)x^*(t)\}$$where $\Bbb E(\cdot)$ denotes the mathematical expectation.

If x(t) is a deterministic signal, then its autocorrelation function will be $$R_x(t)=x(t)*x^*(-t)=\int_{-\infty}^\infty x(\tau)x^*(t+\tau)d\tau$$
where $\cdot *\cdot$ denotes the convolution operation.

Based on the previous category, we have four possible cases:


  (1) Both $x(t)$ and $h(t)$ are deterministic.


In this case, $y(t)$ is deterministic and we obtain $$R_y(t){=\int_{-\infty}^\infty y(\tau)y^*(t+\tau)d\tau\\=\int_{-\infty}^\infty \int_{-\infty}^\infty\int_{-\infty}^\infty x(\tau_1)h(\tau-\tau_1)x^*(\tau_2)h^*(t+\tau-\tau_2)d\tau_1d\tau_2d\tau\\=\int_{-\infty}^\infty \int_{-\infty}^\infty x(\tau_1)x^*(\tau_2)\int_{-\infty}^\infty h(\tau-\tau_1)h^*(t+\tau-\tau_2)d\tau d\tau_1d\tau_2\\=\int_{-\infty}^\infty \int_{-\infty}^\infty x(\tau_1)x^*(\tau_2)\int_{-\infty}^\infty h(\tau)h^*(t+\tau+\tau_1-\tau_2)d\tau d\tau_1d\tau_2\\=\int_{-\infty}^\infty \int_{-\infty}^\infty x(\tau_1)x^*(\tau_2)R_h(t+\tau_1-\tau_2) d\tau_1d\tau_2\\=\int_{-\infty}^\infty \int_{-\infty}^\infty x(\tau_1+\tau_2)x^*(\tau_2)R_h(t+\tau_1) d\tau_1d\tau_2\\=\int_{-\infty}^\infty R_h(t+\tau_1) \int_{-\infty}^\infty x(\tau_1+\tau_2)x^*(\tau_2) d\tau_2 d\tau_1\\=\int_{-\infty}^\infty R_h(t+\tau_1) \int_{-\infty}^\infty x(\tau_2)x^*(\tau_2-\tau_1) d\tau_2 d\tau_1\\=\int_{-\infty}^\infty R_h(t+\tau_1) R_x(-\tau_1) d\tau_1\\=R_x(t)*R_h(t)
}$$


  (2) $x(t)$ is stochastic and $h(t)$ is deterministic.


In this case, $y(t)$ is stochastic and we obtain:
$$R_y(t_1,t_2){=\Bbb E\{y(t_1)y^*(t_2)\}
\\=\Bbb E\left\{\int_{-\infty}^\infty\int_{-\infty}^\infty x(\tau_1)h(t_1-\tau_1)x^*(\tau_2)h^*(t_2-\tau_2)d\tau_1d\tau_2\right\}
\\=\int_{-\infty}^\infty\int_{-\infty}^\infty \Bbb E\{x(\tau_1)x^*(\tau_2)\}h(t_1-\tau_1)h^*(t_2-\tau_2)d\tau_1d\tau_2
\\=\int_{-\infty}^\infty\int_{-\infty}^\infty R_x(\tau_1,\tau_2)h(t_1-\tau_1)h^*(t_2-\tau_2)d\tau_1d\tau_2
\\=\int_{-\infty}^\infty\left[\int_{-\infty}^\infty R_x(\tau_1,\tau_2)h(t_1-\tau_1)d\tau_1\right]h^*(t_2-\tau_2)d\tau_2
\\=\int_{-\infty}^\infty\left[ R_x(t_1,\tau_2)*_{t_1}h(t_1)\right]h^*(t-\tau_2)d\tau_2
\\=\left[ R_x(t_1,t_2)*_{t_1}h(t_1)\right]*_{t_2}h^*(t_2)
}$$where $*_u$ means convolution w.r.t. $u$.


  (3) $h(t)$ is stochastic and $x(t)$ is deterministic.


This case is very similar to case (2). Since the system is LTI, we can reverse the impulse response with input signal to obtain $$R_y(t_1,t_2)=\left[ R_h(t_1,t_2)*_{t_1}x(t_1)\right]*_{t_2}x^*(t_2)$$


  (4) Both $x(t)$ and $h(t)$ are stochastic.


As the statistical sources of $x(t)$ and $h(t)$ are often independent, it is quite reasonable to assume that $x(t)$ and $h(t)$ are also independent. To use this fact, we start by:$$R_y(t_1,t_2){=\Bbb E\{y(t_1)y^*(t_2)\}
\\=\Bbb E\left\{\int_{-\infty}^\infty\int_{-\infty}^\infty x(\tau_1)h(t_1-\tau_1)x^*(\tau_2)h^*(t_2-\tau_2)d\tau_1d\tau_2\right\}
\\=\Bbb \int_{-\infty}^\infty\int_{-\infty}^\infty \Bbb E\left\{x(\tau_1)h(t_1-\tau_1)x^*(\tau_2)h^*(t_2-\tau_2)\right\}d\tau_1d\tau_2
\\=\Bbb \int_{-\infty}^\infty\int_{-\infty}^\infty \Bbb E\{x(\tau_1)x^*(\tau_2)\}\Bbb E\{ h(t_1-\tau_1)h^*(t_2-\tau_2)\}d\tau_1d\tau_2
\\=\Bbb \int_{-\infty}^\infty R_x(t_1,\tau_2)*_{t_1}R_h(t_1,t_2-\tau_2)d\tau_2
\\=R_x(t_1,t_2)*_{t_1}*_{t_2}R_h(t_1,t_2)
}$$
which means that a dual convolution must be performed: one on $t_1$ and one on $t_2$.