Can I replace the distribution in Gaussian Process Regression with a different regression? Sorry for the confusing title. Let me try to clarify:
I have a time series of wind speeds with some missing points here and there. I want to interpolate these points and have tried mainly polynomial and spline interpolation (also weighted interpolation), but for some cases, these just did not work well. Moreover, these interpolation methods do not take into account what a reasonable wind speed value range would be.
Looking into Gaussian Process Regression (GPR), that seems to basically do what I want and uses an underlying Gaussian distribution for the regression / interpolation.
What I am now wondering is if it is possible to use the underlying principle of GPR with a distribution different than Gaussian distribution (e.g. Weibull distribution would make a lot of sense for wind speeds). If I understood it correctly, the used kernel is independant of the distribution, so the kernel should not be a problem here.
So, would it be possible to simply replace the distribution for this method? Or is it simply not possible for some mathematical reason? I have a quite weak mathematical understanding of the GPR, so that is entirely possible. The Wiki article on the topic was pretty hard to understand for me, so I would appretiate it if you could include a simple-to-understand explanation in your answers.
Thanks a lot in advance!
 A: Disclaimer: This is an unfinished answer which still needs some (mathematical / research / googling) work. Feel free to downvote if there is somebody with more experience concerning Weibull processes...

So, would it be possible to simply replace the distribution for this method? Or is it simply not possible for some mathematical reason?

Let us quickly review how Gaussian Processes work. The first assumption is that there is some space $T$ from which the indices of random variables $(X_t)_{t \in T}$ come. In your case, $T$ is probably a series of times (like $T = \{0,1,2,...\}$) and $X_t$ resembles the windspeed at day/hour/second $t$.
In Gaussian processes one assumes that for any finite subset $t_1, ..., t_n$ of times, $(X_{t_1}, ..., X_{t_n})$ is distributed with a density being a multivariate normal distribution $\mathcal{N}(\Sigma, \mu)$ with covariance matrix $\Sigma = (k(t_i, t_j)_{i,j=1,...,n})$ and mean $\mu = (\mu(t_i))_{i=1,...,n}$. So far there is nothing special about this multivariate distribution (i.e. up until now we could easily replace it by another multivariate distribution).
What is so special about the multivariate normal distribution / why don't people use another one?
Well, the unique property of the multivariate normal distribution is that given a 'new' $t_{n+1}$, we can find an explicit expression for the density
  $$f_{X_{t_{n+1}}|X_{t_1}, ..., X_{t_{n}}}(x_{t_{n+1}}|x_{t_1}, ..., x_{t_{n}}) = \frac{f_{X_{t_{n+1}},X_{t_1}, ..., X_{t_{n}}}(x_{t_{n+1}},x_{t_1}, ..., x_{t_{n}})}{f_{X_{t_1}, ..., X_{t_{n}}}(x_{t_1}, ..., x_{t_{n}})}$$
because using a somewhat involved computation one can show the following:
If the whole covariance matrix of $(X_{t_{n+1}},X_{t_{1}}, ..., X_{t_{n}})$ matrix looks like
$$\Sigma_{\text{all}} = \begin{pmatrix} \Sigma_{n+1,n+1} & \Sigma_{n+1,*}^T \\ \Sigma_{n+1,*} & \Sigma_{*,*} \end{pmatrix} = \begin{pmatrix} k(t_{n+1}, t_{n+1}) & k(t_{n+1}, *)^T \\ k(t_{n+1}, *) & k(*,*)\end{pmatrix}$$
i.e.
$$\Sigma_{n+1,*}^T = (k(t_{n+1}, t_1), ..., k(t_{n+1}, t_n))$$
(the covariance pairing between the 'new' point $x_{t_{n+1}}$ and the 'old'/'known' ones $x_{t_1}, ..., x_{t_n}$)
and
$$\Sigma_{*,*} = (k(t_i, t_j))_{i,j=1,...,n}$$
(the covariance matrix on the 'old'/'known' points $x_{t_1}, ..., x_{t_n}$
and if we write the mean of all variables as
$$\mu = (\mu(t_{n+1}), \mu_*)$$
where $\mu_* = (\mu(t_1), ..., \mu(t_n))$ 
then
$f_{X_{t_{n+1}}|X_{t_1}, ..., X_{t_{n}}}$ is a multivariate normal distribution again and 
$$f_{X_{t_{n+1}}|X_{t_1}, ..., X_{t_{n}}} = \mathcal{N}(\Sigma', \mu')$$
with
$$\Sigma' = k(t_{n+1}, t_{n+1}) - \Sigma_{n+1,*}^T (\Sigma_{*,*})^{-1} \Sigma_{n+1,*}$$
and
$$\mu' = \mu(t_{n+1}) + \Sigma_{n+1,*}^T (\Sigma_{*,*})^{-1}((x_{t_1}, ..., x_{t_n})^T - \mu_*)$$
In one sentence: 
What is so special about the multivariate normal distribution / why don't people use another one?
The predictive distribution still belongs to the same class of distributions and given the values $x_{t_1}, ..., x_{t_n}$ and the new point $x_{t_{n+1}}$ we can explicitly compute the parameters (and therefor the distribution itself).
Things you have to do if you want to do this with the Weibull distribution:


*

*First of all make sure that there is a multivariate Weibull distribution. this seems to indicate that there is such a thing and how to compute it.

*Carry out the involved computation to verify that f_{X_{t_{n+1}}|X_{t_1}, ..., X_{t_{n}}} is a univariate Weibull distribution with parameters that can be explicitly computed from the values $x_{t_1}, ..., x_{t_n}, x_{t_{n+1}}$ and the parameters of the distribution of $x_{t_1}, ..., x_{t_n}$.


However, this might be math heavy (I did not try it yet) and it might be the case that somebody already did this...
Hope this helps...
PS: Even I do not understand what this Kriging is about and how exactly it is related to Gaussian Processes from the Wikipedia page :-)
