2

我有一个空间二维域,比如 [0,1]×[0,1]。在这个域中,有 6 个点可以观察到一些感兴趣的标量(例如,温度、机械应力、流体密度等)。如何预测未观察点的兴趣数量?换句话说,我如何在 Python 中插入空间数据?

例如,考虑 2D 域(输入)中的点的以下坐标和感兴趣数量的相应观察值(输出)。

import numpy as np
coordinates = np.array([[0.0,0.0],[0.5,0.0],[1.0,0.0],[0.0,1.0],[0.5,1.],[1.0,1.0]])
observations = np.array([1.0,0.5,0.75,-1.0,0.0,1.0])

X 和 Y 坐标可以通过以下方式提取:

x = coordinates[:,0]
y = coordinates[:,1]

以下脚本创建一个散点图,其中黄色(分别为蓝色)代表高(分别为低)输出值。

import matplotlib.pyplot as plt
fig = plt.figure()
plt.scatter(x, y, c=observations, cmap='viridis')
plt.colorbar()
plt.show()

6 点的观察

我想使用克里金法预测二维输入域内规则网格上的感兴趣标量。我怎样才能在 Python 中做到这一点?

4

1 回答 1

5

OpenTURNS中,KrigingAlgorithm该类可以根据特定输入点的已知输出值估计高斯过程模型的超参数。然后,的getMetamodel方法KrigingAlgorithm返回一个对数据进行插值的函数。

首先,我们需要将 Numpy 数组转换coordinatesobservationsOpenTURNSSample对象:

import openturns as ot
input_train = ot.Sample(coordinates)
output_train = ot.Sample(observations, 1)

数组coordinates有 shape (6, 2),所以它变成了Sample大小为 6 和维度 2 的 a。数组observations有 shape (6,),这是不明确的:它是Sample大小为 6 和维度 1 的 a,还是Sample大小为 1 和维度 6 的 a?Sample为了澄清这一点,我们在对类构造函数的调用中指定维度 (1) 。

下面,我们定义一个具有恒定趋势函数和平方指数协方差核的高斯过程模型:

inputDimension = 2
basis = ot.ConstantBasisFactory(inputDimension).build()
covariance_kernel = ot.SquaredExponential([1.0]*inputDimension, [1.0])
algo = ot.KrigingAlgorithm(input_train, output_train,
                           covariance_kernel, basis)

然后我们拟合趋势的值和协方差核的参数(幅度参数和尺度参数),得到一个元模型:

# Fit
algo.run()
result = algo.getResult()
krigingMetamodel = result.getMetaModel()

结果krigingMetamodel是 a Function,它以 2DPoint作为输入并返回 1D Point。它预测感兴趣的数量。为了说明这一点,让我们构建二维域 [0,1]×[0,1] 并用规则网格对其进行离散化:

# Create the 2D domain
myInterval = ot.Interval([0.0, 0.0], [1.0, 1.0])
# Define the number of interval in each direction of the box
nx = 20
ny = 10
myIndices = [nx - 1, ny - 1]
myMesher = ot.IntervalMesher(myIndices)
myMeshBox = myMesher.build(myInterval)

使用我们krigingMetamodel来预测该网格上感兴趣的数量所取的值可以通过以下语句完成。我们首先vertices将网格的 设为Sample,然后通过对元模型的一次调用来评估预测(这里不需要for循环):

# Predict
vertices = myMeshBox.getVertices()
predictions = krigingMetamodel(vertices)

为了使用 Matplotlib 查看结果,我们首先必须创建pcolor函数所需的数据:

# Format for plot
X = np.array(vertices[:, 0]).reshape((ny, nx))
Y = np.array(vertices[:, 1]).reshape((ny, nx))
predictions_array = np.array(predictions).reshape((ny,nx))

以下脚本生成绘图:

# Plot
import matplotlib.pyplot as plt
fig = plt.figure()
plt.pcolor(X, Y, predictions_array)
plt.colorbar()
plt.show()

克里金元模型的预测

我们看到元模型的预测等于观察到的输入点的观察。

这个元模型是坐标的平滑函数:它的平滑度随着协方差核平滑度的增加而增加,平方指数协方差核恰好是平滑的。

于 2019-12-13T22:19:03.907 回答