GaussianProcessRegressor: wrong std and cov results when n_features>1 and no y normalization

### Describe the bug

When `n_features > 1` and `normalization_y` is `False`, the `GaussianProcessRegressor.predict` seems to return bad std and cov results, as it doesn't consider the scale of the different features (while it seems to be ok when `n_features > 1` and `normalization_y` is `True`).

By taking a look at the code, we can see that `GaussianProcessRegressor.predict` uses the `_y_train_std` attribute to compute the variance and covariance but this attribute is set to `ones(n_features)` when `normalize_y` is set to `False` (default value), giving equal scale to all features.

To fix this bug, one should always compute `_y_train_std` from the training data and use the boolean attribute `normalize_y` to undo the normalization of `y_mean` if necessary.

### Steps/Code to Reproduce

```python
import pytest
from numpy import array
from numpy import hstack
from sklearn.gaussian_process import GaussianProcessRegressor

x = array([[0.], [0.5], [1.]])
y = hstack((x**2, 10*x**2))
# Note that the second output is equal to 10 times the first one.

# With output normalization
gpr = GaussianProcessRegressor(normalize_y=True)
gpr.fit(x, y)
std = gpr.predict(array([[0.25]]), return_std=True)[1][0]
assert std[0] != std[1]
assert std[0] == pytest.approx(std[1]/10, rel=1e-9)
# As expected, the variance of the second output is 10 times larger than the first output.

# Without output normalization
gpr = GaussianProcessRegressor(normalize_y=False)
gpr.fit(x, y)
std = gpr.predict(array([[0.25]]), return_std=True)[1][0]
assert std[0] == std[1]
# The variance of the second output is equal to the variance of the first output.
```

### Expected Results

Without output normalization, the variance of the second output should be 10 times larger than the first output.

### Actual Results

Without output normalization, the variance of the second output is equal to the variance of the first output.

### Versions

```shell
System:
    python: 3.9.12 (tags/v3.9.12:b28265d, Mar 23 2022, 23:52:46) [MSC v.1929 64 bit (AMD64)]
executable: C:\Users\matthias.delozzo\workspace\GEMSEO\gemseo\.tox\py39\.venv\Scripts\python.exe
   machine: Windows-10-10.0.19045-SP0

Python dependencies:
      sklearn: 1.5.1
          pip: None
   setuptools: None
        numpy: 1.26.4
        scipy: 1.13.1
       Cython: None
       pandas: 2.2.2
   matplotlib: 3.9.2
       joblib: 1.4.2
threadpoolctl: 3.5.0

Built with OpenMP: True

threadpoolctl info:
       user_api: blas
   internal_api: openblas
    num_threads: 8
         prefix: libopenblas
       filepath: C:\Users\matthias.delozzo\workspace\GEMSEO\gemseo\.tox\py39\.venv\Lib\site-packages\numpy.libs\libopenblas64__v0.3.23-293-gc2f4bdbb-gcc_10_3_0-2bde3a66a51006b2b53eb373ff767a3f.dll
        version: 0.3.23.dev
threading_layer: pthreads
   architecture: Haswell

       user_api: openmp
   internal_api: openmp
    num_threads: 8
         prefix: vcomp
       filepath: C:\Users\matthias.delozzo\workspace\GEMSEO\gemseo\.tox\py39\.venv\Lib\site-packages\sklearn\.libs\vcomp140.dll
        version: None

       user_api: blas
   internal_api: openblas
    num_threads: 8
         prefix: libopenblas
       filepath: C:\Users\matthias.delozzo\workspace\GEMSEO\gemseo\.tox\py39\.venv\Lib\site-packages\scipy.libs\libopenblas_v0.3.27--3aa239bc726cfb0bd8e5330d8d4c15c6.dll
        version: 0.3.27
threading_layer: pthreads
   architecture: Haswell
```

COMMENT:
The code I posted fixes the reported bug, running the code to reproduce the error we get that _y_train_std reflects the true scale of the features, even when normalization is disabled. 

However, this change causes some test failures because the predict function in GaussianProcessRegressor currently assumes normalization is always applied when undoing normalization during prediction.

I belive that, to fully resolve the issue, the predict function should be modified to conditionally apply the normalization logic based on the normalize_y flag.

I will see if in the next days I manage to do that.

Labels: Bug