VGP stddev not working properly

Hallo,

I'm using the variational gaussian process fitting from tf.probability and so fare everything works fine. The mean() on the final vgp gives me for predictions a nice function but when I try to plot also the stddev, this works not with all my datasets. Below are two examples, where I plotted the means and the stddev. In the second one one can see, that it plots the stddev, but with some spikes and I think the plotted stddev is too large (I think the spikes are at the location of the inducing points). Such spikes are also visible in the mean. But this spikes are not visible at locations, where no training-data were located, so when I realy make a forecast for example. But in the first example there is never even a stddev visible (or it is very small).

Is it a bug or did I handle something wrong?

Following my definition of the Kernel:

amplitude = tf.exp(tf.Variable(1, dtype=dtype, name='amplitude', use_resource=True))
length_scale = tf.exp(tf.Variable(3, dtype=dtype, name='length_scale', use_resource=True))
smother_kernel = tfk.ExponentiatedQuadratic(amplitude=amplitude, length_scale=length_scale)

# periodic kernel hyperparameters:
periodic_amplitude_1y = tf.exp(tf.Variable(5, dtype=dtype, name='periodic_amplitude_1y', use_resource=True))
periodic_length_scale_1y = tf.exp(tf.Variable(0.69, dtype=dtype, name='periodic_length_scale_1y', use_resource=True))
periodic_period_1y = tf.exp(tf.Variable(9.078, dtype=dtype, name='periodic_period_1y', use_resource=True))  # period of 24h * 365d

# periodic kernel:
local_periodic_kernel = (tfk.ExpSinSquared(amplitude=periodic_amplitude_1y,
                                           length_scale=periodic_length_scale_1y,
                                           period=periodic_period_1y)
                         )

# linear kernel parameters:
linear_bias_variance = tf.Variable(0.00001, dtype=dtype, name='linear_bias_variance')  # needs a small start value, otherwise it is not trained (if too large, it crashes!)
linear_slope_variance = tf.Variable(0.00001, dtype=dtype, name='linear_slope_variance')  # needs a small start value, otherwise it is not trained (if too large, it crashes!)
linear_shift = tf.Variable(0.00001, dtype=dtype, name='linear_slope_variance')  # needs a small start value, otherwise it is not trained (if too large, it crashes!)

# linear kernel for trends (pollution of solar-modules):
linear_kernel = tfk.Linear(
    bias_variance=linear_bias_variance,
    slope_variance=linear_slope_variance,
    shift=linear_shift
)

kernel = (linear_kernel + local_periodic_kernel + smother_kernel)

observation_noise_variance = tf.exp(tf.Variable(0.5, dtype=dtype, name='observation_noise_variance', use_resource=True))

No source file or test is named. Start by reproducing the reported behavior with tf.probability's variational Gaussian process, the supplied kernel definition, and both datasets, then inspect how predictions and standard deviations are computed around inducing points. Done means determining whether the spikes or near-zero standard deviation are expected or a defect, with a confirmed fix or clear usage explanation.