Find optimal parameters for UK with external drifts

[jlemond]

Hello,

I currently use your very useful kriging package PyKrige for Universal Kriging with 2 external drifts (gridded temperature and gridded elevation) as follow:

UK = UniversalKriging(
    lon,
    lat,
    temp,
    variogram_model=variogram,
    nlags=nlags,
    exact_values=True,
    verbose=True,
    enable_plotting=True,
    drift_terms=['specified'],
    specified_drift=[z,raw],
)
 
krigresult,var =UK.execute('grid',gridx,gridy,specified_drift_arrays=[dem_indo,raw_temp_fit.values])

where:
lon= longitudes of stations,
lat= latitudes of stations,
temp= the temperature values of stations that I want to krige,
z= the elevation of stations used as external drift,
raw= the temperature values used as external drift,
gridx= vector of longitudes,
gridy= vector of latitudes ,
dem_indo= gridded elevation data ,
raw_temp_fit.values= gridding temperature data

I would like to use Kriging CV to search optimal parameters amongst different variograms and number of bins (nlags).

I’ve tried following code:

param_dict= {
    "method": ['universal'],
    "variogram_model": ["linear", "power", "gaussian", "spherical", "exponential"],
    "nlags": [6, 10, 15, 20,30],
    "drift_terms":['specified']
}
    
estimator = GridSearchCV(Krige(), param_dict, verbose=True, return_train_score=True)
 
X=np.array([lon,lat]); X=X.T
y=val    
estimator.fit(X=X, y=val)

As I don’t find the specified drift argument in the Krige() function, I wonder how to manage the drift values. Could you help me?

Thanks in advance,
Julien

[MuellerSeb]

The Krige class doesn't support external drift kriging at the moment.

I extended the original class last year to be able to use universal kriging, but since external drift kriging needs a specified drift at the output locations, I didn't know how to properly pass them.

This was decided here: #158 (comment)

We need further discussion about that internally.

[MuellerSeb]

For now you could use ext_drift_grid of the Krige class And use a grid defined with the external drift, where all values are extracted during calculation. But that only works for 2D.

See:
https://geostat-framework.readthedocs.io/projects/pykrige/en/stable/generated/pykrige.rk.Krige.html#pykrige.rk.Krige