如何在 openturns 中获得更好的克里金结果图?

Posted

技术标签:

【中文标题】如何在 openturns 中获得更好的克里金结果图?【英文标题】:how to get better Kriging result graphs in openturns? 【发布时间】:2021-10-25 18:01:31 【问题描述】:

我执行了球面克里金法,但我似乎无法获得好的输出图。 坐标(x 和 y)的范围从大约 51 纬度到大约 6.5 经度 我的观察范围从 -70 到 +10 这是我的代码:

import openturns as ot
import pandas as pd
# your input / output data can be easily formatted as samples for openturns


df = pd.read_csv("kreuzkerpenutm.csv")


inputdata = ot.Sample(df[['x','y']].values)
outputdata = ot.Sample(df[['z']].values)


dimension = 2  # dimension of your input (x,y)
basis = ot.ConstantBasisFactory(dimension).build()
covarianceModel = ot.SphericalModel(dimension)
    
algo = ot.KrigingAlgorithm(inputdata, outputdata, covarianceModel, basis)
algo.run()
result = algo.getResult()
metamodel = result.getMetaModel()

lower = [-10.0] * 2 # lower bound of the 2D window
upper = [50.0] * 2 # upper bound of the 2D window
graph = metamodel.draw(lower, upper)
graph.setBoundingBox(ot.Interval(lower, upper))
graph.add(ot.Cloud(inputdata)) # overlay a scatter plot of the observation points
graph.setTitle("Kriging metamodel")

# A View object allows us to interact with the underlying matplotlib figure 
from openturns.viewer import View
view = View(graph, legend_kw='bbox_to_anchor':(1,1), 'loc':"upper left")
view.getFigure().tight_layout()

这是我的输出:

kriging metamodel graph

我不知道为什么我的图表不会显示我的输入以及我的克里金结果。

感谢您的想法和帮助

【问题讨论】:

提供和回答并不容易,因为我无法访问 .csv 文件。您会将脚本更新为重现您的问题的最小工作示例吗? 【参考方案1】:

如果输入数据未按 [-1,1]^d 进行缩放,则克里金元模型可能无法使用最大似然优化来识别缩放参数。为了对此提供帮助,我们可以:

为协方差模型的尺度参数提供更好的起点(这是下面的技巧“A”), 设置优化算法的界限,以便搜索参数的区间与手头的数据相对应(这是下面的技巧“B”)。

这就是以下脚本的作用,使用模拟数据而不是 csv 数据文件。在脚本中,我使用 g 函数创建数据,该函数被缩放,以便它产生 [-10, 70] 范围内的结果,就像您的问题一样。请仔细查看设置协方差模型初始值的setScale() 方法:这是优化算法的起点。然后看setOptimizationBounds()方法,它设置了优化算法的界限。

import openturns as ot


dimension = 2  # dimension of your input (x,y)
distribution = ot.ComposedDistribution([ot.Uniform(-10.0, 50.0)] * dimension)
inputdata = distribution.getSample(100)

g = ot.SymbolicFunction(["x", "y"], ["30 + 3.0 * sin(x / 10.0) * (y / 10.0) ^ 2"])

outputdata = g(inputdata)

basis = ot.ConstantBasisFactory(dimension).build()
covarianceModel = ot.SphericalModel(dimension)
covarianceModel.setScale(inputdata.getMax())  # Trick A
algo = ot.KrigingAlgorithm(inputdata, outputdata, covarianceModel, basis)
# Trick B, v2
x_range = inputdata.getMax() - inputdata.getMin()
scale_max_factor = 2.0  # Must be > 1, tune this to match your problem
scale_min_factor = 0.1  # Must be < 1, tune this to match your problem
maximum_scale_bounds = scale_max_factor * x_range
minimum_scale_bounds = scale_min_factor * x_range
scaleOptimizationBounds = ot.Interval(minimum_scale_bounds, maximum_scale_bounds)
algo.setOptimizationBounds(scaleOptimizationBounds)
algo.run()
result = algo.getResult()
metamodel = result.getMetaModel()
metamodel.setInputDescription(["x", "y"])
metamodel.setOutputDescription(["z"])

lower = [-10.0] * 2  # lower bound of the 2D window
upper = [50.0] * 2  # upper bound of the 2D window
graph = metamodel.draw(lower, upper)
graph.setBoundingBox(ot.Interval(lower, upper))
graph.add(ot.Cloud(inputdata))  # overlay a scatter plot of the observation points
graph.setTitle("Kriging metamodel")

# A View object allows us to interact with the underlying matplotlib figure
from openturns.viewer import View

view = View(graph, legend_kw="bbox_to_anchor": (1, 1), "loc": "upper left")
view.getFigure().tight_layout()

前面的脚本生成下图。

还有其他方法可以实现技巧 B。这是 J.Pelamatti 提供的一种方法:

# Trick B, v3
for d in range(X_train.getDimension()):
    dist = scipy.spatial.distance.pdist(X_train[:,d])
    scale_max_factor = 2.0  # Must be > 1, tune this to match your problem
    scale_min_factor = 0.1  # Must be < 1, tune this to match your problem
    maximum_scale_bounds = scale_max_factor * np.max(dist)
    minimum_scale_bounds = scale_min_factor * np.min(dist)

这个话题在this particular thread in OT's forum讨论。

【讨论】:

【参考方案2】:

抱歉回复晚了。

您使用的是哪个版本的openturns? 可能您对(输入)数据进行了嵌入式转换,这使得数据范围大致介于 (-3, 3) 之间(标准缩放)。在这种情况下,克里金结果应该包含转换。 随着最近的openturns implementations,此功能已被删除。

希望这会有所帮助。 干杯

【讨论】:

这并没有提供问题的答案。一旦你有足够的reputation,你就可以comment on any post;相反,provide answers that don't require clarification from the asker。 - From Review

以上是关于如何在 openturns 中获得更好的克里金结果图?的主要内容,如果未能解决你的问题,请参考以下文章

数据可视化应用Python-pykrige包-克里金(Kriging)插值计算及可视化绘制

如何在 Python 中使用克里金法插入测站数据?

如何在 Python 中使用克里金法对 2D 空间数据进行插值?

如何在大型数据集上执行克里金(高斯过程回归)?

自动克里金法和 SpatialPixel 网格文件

R: Kriging interpolation and cross validation 克里金插值及交叉验证浅析