The sampling offset is a time delay between the I and Q paths. The Fourier Transform of a time delay $T$ is $e^{-j\omega T}$ and thus we see that this will impart a linear phase that is proportional to frequency (and knowing $T$ and our frequencies of interest, we know exactly what this phase is). In contrast to the Hilbert Transform, which has a constant quadrature phase across all frequencies, the phase due to the time offset is varying linearly over frequency. If the sampling rate is high enough, and the signal bandwidth and spectral occupancy low enough, the phase may be insignificant and can be ignored, but in most practical cases where we are not significantly oversampling the signal, this added delay must be compensated for. This is easily done with digital techniques such as Pade approximations, farrow filters, polyphase filters etc.
Question 1:
In this case of the 48 KHz sampling with every other sample selecting either I or Q from the quadrature detector, we note that the time offset between the I and Q paths with each sampled at 24 KHz is 1/48 KHz, and given the Fourier relationship above, this would impart a phase error from true quadrature proportional to frequency as given by $2 \pi f/48e3$. The 4 KHz component of the waveform would therefore have a phase error of $2\pi 4/48 = 0.52 $ radians. This would lead to a significant sideband error (0.52 radians is no longer small angle, but for small angles the relationship between quadrature error and sideband suppression is $20log_{10}(\phi/2)$.) A chart showing the sideband level versus phase and amplitude error is available here, where we see that 0.52 radians (or 29.8 degrees) is far off the chart; meaning unacceptable quadrature error.
Question 2:
Zero stuffing replicates the spectrum without introducing any other distortion within the bandwidth of the original spectrum, but adds spectral distortion out of the bandwidth through spectral images. For an upsampling by two of a 24 KHz waveform with a 4 KHz signal bandwidth - we would get the following spectrum from the zero-stuffing prior to subsequent interpolation filtering:

The upper spectrum is the original spectrum for each I or Q path sampled at 24 KHz, out to the Nyquist frequency of $f_s/2 = 24$ KHz. The lower spectrum shows the spectral replication that occurs due to the zero-stuffing: the Nyquist frequency extends out to twice the original or 24 KHz, but now an image of the original spectrum exists right near the Nyquist boundary. The "perfect" interpolation filter would pass the desired spectrum out to 4KHz with no distortion and completely eliminate the image. Such "perfect" filters are not possible to implement, but we can constrain the filter complexity based on the allowable distortion we can accept.
Using linear interpolation to replace the zeros is identical to using an FIR filter with coefficients given as [0.5, 1, 0.5]. The frequency response for this filter is shown below:

This filter adds about 0.7 dB of passband droop at 4 KHz, and attenuates the image starting at 20 KHz by 24 dB, with increasing attenuation after that.
In comparison below I show what can be achieved with a simple halfband FIR filter using the least-squares algorithm for minimum error in the 0 to 4 KHz passband and maximum rejection in the 20 KHz to 24 KHz rejection band. I designed the filter for 7 taps, but as we'll see further below this ends up being an interpolator using 2 pairs of adjacent samples instead of 1 pair as the OP is doing with the linear interpolator. Here we see for comparison to the linear interpolator how we can get greater than 50 dB rejection at the starting boundary of the image spectrum for this case.

The normalized coefficients for this filter are:
h[0] = h[6] = -0.07454
h[1] = h[5] = 0
h[2] = h[4] = 0.57329
h[3] = 1
So for a 7 sample buffer with samples numbered 0 through 6; add samples 0 and 6 together and multiply by -0.07454, add that to the sum of samples 2 and 4 multiplied by 0.57329, and add that to sample 3. Since the buffer will have alternating zeros from the stuffing, what ends up happening in the operation as a filter is that when the filter is centered on a non-zero value, all the other coefficients will be zero'd and only that one value is selected (as we would want), and when the filter is centered on a zero'd value (as in the sketch below) the interpolated sample is found from just two sums and two products of the two adjacent non-zero samples on each side.

I used the notation $x_{-1}$ to indicate the first earlier sample and $x_{+1}$ for the first later sample, relative to the zero being replaced with an interpolated value.
This quantifies what that distortion would be for the OP's case between a simple linear interpolation versus a more sophisticated filter, and I have provided an example filter that could alternatively be used with significantly better rejection of the image distortion and significantly less passband droop, with the added complexity of using 2 adjacent samples on each side instead of one. Which way to go depends on the requirements for what level of distortion would be acceptable.
It was asked in the comments how a CIC would compare so I also included the frequency responses for a 1 section CIC (or first order CIC) as "CIC", and cascaded CICs as "CIC^2", "CIC^3", and "CIC^5". A 2 section (second order) CIC would have the same frequency response as the linear interpolator the OP suggested. Order is increased by simply including additional "combs" and "integrator" stages as diagramed below for an Nth order CIC (CIC^N):

In order to get greater than 50 dB rejection (for comparison to the improved interpolator I provided), a 5 section (fifth order) CIC would be needed. I didn't show the fourth order case; this one would have -47 dB rejection at the 20 KHz edge of the image location. But also note the significant increase in passband droop from DC to 4 KHz that the CIC filter would give. The passband droop if an issue can easily be compensated with an additional 3 tap FIR filter, but at this point the solution isn't nearly as simple compared to the enhanced interpolator based on a halfband FIR filter, especially if you are working with a fixed point implementation as a higher order CIC will easily introduce overflow errors if precision isn't carefully managed.
