求助：使用Kriging插值绘制阿根廷冰雹报告地图（R/Python）

我尝试用Kriging插值法绘制阿根廷冰雹报告的分布地图，但在R语言中无法成功实现，希望能得到技术支持。

我的冰雹报告数据包含日期、经度（Longitud）、纬度（Latitud）、冰雹发生记录标识等字段，为点格式的观测记录。

我借助ChatGPT编写了如下R代码，但运行后无法得到预期的插值地图：

library(gstat)
library(sp)

# 加载包含冰雹报告数据的CSV文件
datos <- read.csv("C:/Users/Usuario/Desktop/Agus/FACULTAD/Tesis/SI - NO ROMINA.csv", nrows = 7089)  # 替换为你的CSV文件路径及名称

# 将纬度和经度列转换为数值格式
datos$Latitud <- as.numeric(as.character(datos$Latitud))
datos$Longitud <- as.numeric(as.character(datos$Longitud))

# 创建包含经纬度数据的SpatialPoints对象
datos_sp <- SpatialPointsDataFrame(
coords = datos[, c("Longitud", "Latitud")],
data = datos,
proj4string = CRS("+proj=longlat +datum=WGS84")
)

# 创建用于预测的规则网格
grd <- expand.grid(
Longitud = seq(min(datos$Longitud, na.rm = TRUE), max(datos$Longitud, na.rm = TRUE), length.out = 100),
Latitud = seq(min(datos$Latitud, na.rm = TRUE), max(datos$Latitud, na.rm = TRUE), length.out = 100)
)
grd_sp <- SpatialPoints(coords = grd, proj4string = CRS("+proj=longlat +datum=WGS84"))

# 执行Kriging插值
kriged <- krige(rep(1, nrow(datos_sp)) ~ 1, datos_sp, grd_sp)

# 获取网格各点的预测值
predicted_values <- vector("numeric", length = nrow(grd))
for (i in 1:nrow(grd)) {
predicted_values[i] <- krige(rep(1, nrow(datos_sp)) ~ 1, datos_sp, grd_sp[i, ])$predictions
}

# 将预测值转换为矩阵
kriged_values <- matrix(predicted_values, nrow = 100)

# 绘制Kriging插值地图
par(mar = c(1, 1, 1, 1)) # 必要时调整边距
image(
x = seq(min(datos$Longitud, na.rm = TRUE), max(datos$Longitud, na.rm = TRUE), length.out = 100),
y = seq(min(datos$Latitud, na.rm = TRUE), max(datos$Latitud, na.rm = TRUE), length.out = 100),
z = kriged_values,
col = terrain.colors(50),
main = "Kriging插值地图"
)

请帮忙排查代码问题，或提供可行的Kriging插值实现方案，以生成阿根廷冰雹报告的分布地图。

内容的提问来源于stack exchange，提问作者Agustin Granato