Difference between convolving before/after discretizing LTI systems
Assuming you refer to the ZOH-discretization as shown in this figure from the mathworks site



the relationship between the continuous-time signals and the discrete-time signals can be derived as follows. The signal $u(t)$ is given by

$$u(t)=\sum_ku[k]g(t-kT)\tag{1}$$

where $g(t)$ is a rectangular impulse response (constant in the interval $t\in [0,T]$, zero anywhere else) representing the ZOH, and $T$ is the sampling period. Let $h(t)$ be the impulse response corresponding to the transfer function $H(s)$. Furthermore, let $f(t)$ be the convolution of $g(t)$ and $h(t)$:

$$f(t)=(g\star h)(t)\tag{2}$$

Then the continuous-time output signal $y(t)$ is given by

$$\begin{align}y(t)&amp;=(u\star h)(t)\\
&amp;=\int_{-\infty}^{\infty}u(\tau)h(t-\tau)d\tau\\&amp;=\sum_ku[k]\int_{-\infty}^{\infty}g(\tau-kT)h(t-\tau)d\tau\\&amp;=\sum_ku[k]f(t-kT)\tag{3}\end{align}$$

where I've used Equations $(1)$ and $(2)$. From $(3)$ the discrete-time output signal is easily obtained as

$$y[n]=y(nT)=\sum_ku[k]f((n-k)T)=\sum_ku[k]f[n-k]\tag{4}$$

with $f[n]=f(nT)$. So $y[n]$ is simply the discrete-time convolution of the input signal $u[n]$ and the total impulse response $f[n]$.

From this result it immediately follows that linearity must be satisfied, simply because convolution is a linear operation:

$$f(t)=\alpha_1(g\star h_1)(t)+\alpha_2(g\star h_2)(t)=(g\star (\alpha_1h_1+\alpha_2h_2))(t)\tag{5}$$

However, the (ZOH-)discretization of the concatenation of two transfer functions is not equivalent to the concatenation of the two discretized transfer functions. Referring to the figure above, the first case is equivalent to replacing $H(s)$ by the concatenation of $H_1(s)$ and $H_2(s)$. The other case involves concatenating two complete systems as shown in the figure, the first with $H(s)=H_1(s)$, the second with $H(s)=H_2(s)$. The difference between the two cases is that in the first case you only have one ZOH, whereas in the second one you get two ZOHs. The equivalent discrete-time impulse responses are

$$f[n]=(g\star h_1\star h_2)(nT)\tag{6}$$

and

$$f[n]=(g\star h_1)\star (g\star h_2)(nT)\tag{7}$$

which are generally not identical. The first case in Eq. $(6)$ is a ZOH discretization of the total transfer function $H(s)=H_1(s)H_2(s)$, whereas the second case in Eq. $(7)$ corresponds to a first-order hold (FOH) discretization of the total transfer function.