Anisotropic

In [1]:

from debiased_spatial_whittle.backend import BackendManager

from debiased_spatial_whittle.backend import BackendManager

In [2]:

BackendManager.set_backend("numpy")

BackendManager.set_backend("numpy")

In [3]:

xp = BackendManager.get_backend()

xp = BackendManager.get_backend()

In [4]:

from debiased_spatial_whittle.grids.base import RectangularGrid
from debiased_spatial_whittle.models.univariate import (
    SquaredExponentialModel,
    NuggetModel,
    AnisotropicModel,
)
from debiased_spatial_whittle.sampling.simulation import SamplerOnRectangularGrid
from debiased_spatial_whittle.inference.periodogram import Periodogram, ExpectedPeriodogram
from debiased_spatial_whittle.inference.least_squares import LeastSquareEstimator
from debiased_spatial_whittle.inference.likelihood import DebiasedWhittle, Estimator

from debiased_spatial_whittle.grids.base import RectangularGrid from debiased_spatial_whittle.models.univariate import ( SquaredExponentialModel, NuggetModel, AnisotropicModel, ) from debiased_spatial_whittle.sampling.simulation import SamplerOnRectangularGrid from debiased_spatial_whittle.inference.periodogram import Periodogram, ExpectedPeriodogram from debiased_spatial_whittle.inference.least_squares import LeastSquareEstimator from debiased_spatial_whittle.inference.likelihood import DebiasedWhittle, Estimator

In [5]:

import matplotlib.pyplot as plt

import matplotlib.pyplot as plt

Set up grid and model¶

In [6]:

grid = RectangularGrid((256, 256))
model = SquaredExponentialModel(rho=12.0)
model = AnisotropicModel(model, eta=1.7, phi=xp.pi / 4)
model = NuggetModel(model, nugget=1e-2)
model

grid = RectangularGrid((256, 256)) model = SquaredExponentialModel(rho=12.0) model = AnisotropicModel(model, eta=1.7, phi=xp.pi / 4) model = NuggetModel(model, nugget=1e-2) model

Out[6]:

NuggetModel()

Name	Value	Type	Range
nugget	0.01	ModelParameter	(0, 1)

AnisotropicModel()

Name	Value	Type	Range
eta	1.7	ModelParameter	(0, inf)
phi	0.7853981633974483	ModelParameter	(-1.5707963267948966, 1.5707963267948966)

SquaredExponentialModel()

Name	Value	Type	Range
rho	12.0	ModelParameter	(0, inf)
sigma	1.0	ModelParameter	(0, inf)

Sample from the model¶

In [7]:

sampler = SamplerOnRectangularGrid(model, grid)
data = sampler()

sampler = SamplerOnRectangularGrid(model, grid) data = sampler()

In [8]:

plt.figure()
plt.pcolor(data, cmap="RdBu")
plt.show()

plt.figure() plt.pcolor(data, cmap="RdBu") plt.show()

Set up estimation¶

In [9]:

periodogram = Periodogram()
expected_periodogram = ExpectedPeriodogram(grid, periodogram)

periodogram = Periodogram() expected_periodogram = ExpectedPeriodogram(grid, periodogram)

In [10]:

model_est = SquaredExponentialModel(rho=10.0)
model_est.param.rho.bounds = (5, 100)
model_est.param.sigma.bounds = (0.1, 10)
model_est = AnisotropicModel(model_est)

model_est = SquaredExponentialModel(rho=10.0) model_est.param.rho.bounds = (5, 100) model_est.param.sigma.bounds = (0.1, 10) model_est = AnisotropicModel(model_est)

In [11]:

model_est_ = NuggetModel(model_est, nugget=1e-3)

model_est_ = NuggetModel(model_est, nugget=1e-3)

In [12]:

debiased_whittle = DebiasedWhittle(periodogram, expected_periodogram)

debiased_whittle = DebiasedWhittle(periodogram, expected_periodogram)

Carry out least square fit¶

In [13]:

least_square = LeastSquareEstimator(periodogram, expected_periodogram)
least_square(data, model_est_)
print(model_est_.free_parameter_values_to_array_deep())
model_est_

least_square = LeastSquareEstimator(periodogram, expected_periodogram) least_square(data, model_est_) print(model_est_.free_parameter_values_to_array_deep()) model_est_

[ 9.70342000e-03  5.88931997e-01 -7.76958882e-01  1.21815401e+01
  1.43832414e+00]

Out[13]:

NuggetModel()

Name	Value	Type	Range
nugget	0.009703419998102576	ModelParameter	(0, 1)

AnisotropicModel()

Name	Value	Type	Range
eta	0.5889319967848157	ModelParameter	(0, inf)
phi	-0.7769588817896067	ModelParameter	(-1.5707963267948966, 1.5707963267948966)

SquaredExponentialModel()

Name	Value	Type	Range
rho	12.18154014422066	ModelParameter	(5, 100)
sigma	1.4383241438327155	ModelParameter	(0.1, 10)

Carry out Debiased Whittle fit¶

In [14]:

estimator = Estimator(debiased_whittle)
estimator(model_est_, data)
print(model_est_.free_parameter_values_to_array_deep())
model_est_

estimator = Estimator(debiased_whittle) estimator(model_est_, data) print(model_est_.free_parameter_values_to_array_deep()) model_est_

[ 9.30003207e-03  5.92912190e-01 -7.89126893e-01  1.20333727e+01
  1.03971765e+00]

Out[14]:

NuggetModel()

Name	Value	Type	Range
nugget	0.009300032072528996	ModelParameter	(0, 1)

AnisotropicModel()

Name	Value	Type	Range
eta	0.5929121898050794	ModelParameter	(0, inf)
phi	-0.789126893015422	ModelParameter	(-1.5707963267948966, 1.5707963267948966)

SquaredExponentialModel()

Name	Value	Type	Range
rho	12.0333726920572	ModelParameter	(5, 100)
sigma	1.0397176516824167	ModelParameter	(0.1, 10)