Python进行克里金插值可视化

处理数据并进行克里金插值

# 处理数据后进行克里金插值
import pandas as pd
pm25_region=pd.read_csv(r"C:\Users\LHW\Desktop\Hubei task120210918\extract_pm25\PM25_region.csv")
del pm25_region['type']
year_list=['2019','2020','2021']
data_list=[]
for year in year_list:pm25_region['date'] = pd.to_datetime(pm25_region['date'])pm25_region1=pm25_region.set_index('date')pm2019=pm25_region1[year+'-02':year+'-03']pm2019_mean=pm2019.mean().to_frame()pm2019_mean.columns=[year+'_pm25']data_list.append(pm2019_mean)ini=data_list[0]
ini=ini.join(data_list[1])
ini=ini.join(data_list[2])loc=pd.read_csv("region_aim_loc_info.csv")
loc=loc.set_index('监测点编码')
total_info=loc.join(ini)
del total_info['Unnamed: 0']total_info=total_info.dropna()# 进行克里金插值
import numpy as np
import geopandas as gpd
js=gpd.read_file(r'Hubei.json')
js_box = js.geometry.total_boundsgrid_lon = np.linspace(js_box[0],js_box[2],400)
grid_lat = np.linspace(js_box[1],js_box[3],400)lons = total_info["经度"].values
lats = total_info["纬度"].values
data = total_info["2019_pm25"].valuesfrom pykrige.ok import OrdinaryKriging
OK = OrdinaryKriging(lons, lats, data, variogram_model='gaussian',nlags=6)
z1, ss1 = OK.execute('grid', grid_lon, grid_lat)#转换成网格
xgrid, ygrid = np.meshgrid(grid_lon, grid_lat)
#将插值网格数据整理
df_grid =pd.DataFrame(dict(long=xgrid.flatten(),lat=ygrid.flatten()))
#添加插值结果
df_grid["Krig_gaussian"] = z1.flatten()

站点数据可视化

import plotnine
from plotnine import *
plotnine.options.figure_size = (6, 4)
idw_scatter = (ggplot() +geom_map(js,fill='none',color='gray',size=0.6) +geom_point(total_info,aes(x='经度',y='纬度',fill='2019_pm25'),size=5) +scale_fill_cmap(cmap_name='Spectral_r',name='PM2.5',limits = (0, 100),breaks=[30,40,60,80])+scale_x_continuous(breaks=[108,110,112,114,116,118])+#labs(title="Map Charts in Python Exercise 02: Map IDM point",)+#添加文本信息#annotate('text',x=116.5,y=35.3,label="processed map charts with plotnine",ha="left",size=10)+#annotate('text',x=120,y=30.6,label="Visualization by DataCharm",ha="left",size=9)+theme(text=element_text(family = "SimHei"),#修改背景panel_background=element_blank(),axis_ticks_major_x=element_blank(),axis_ticks_major_y=element_blank(),axis_text=element_text(size=12),axis_title = element_text(size=14,weight="bold"),panel_grid_major_x=element_line(color="gray",size=.5),panel_grid_major_y=element_line(color="gray",size=.5),))
idw_scatter.save('test.jpg',dpi=300)

克里金插值后数据可视化

import plotnine
from plotnine import *
plotnine.options.figure_size = (6.4, 4)
Krig_inter_grid = (ggplot() + geom_tile(df_grid,aes(x='long',y='lat',fill='Krig_gaussian'),size=0.1) +geom_map(js,fill='none',color='gray',size=0.6) + scale_fill_cmap(cmap_name='Spectral_r',name='PM2.5',breaks=[30,40,60,80])+scale_x_continuous(breaks=[108,110,112,114,116,118])+scale_y_continuous(breaks=[28,30,32,34])+#labs(title="Map Charts in Python Exercise 02: Map point interpolation",)+#添加文本信息#annotate('text',x=116.5,y=35.3,label="processed map charts with plotnine",ha="left", size=10)+#annotate('text',x=120,y=30.6,label="Visualization by DataCharm",ha="left",size=9)+theme(text=element_text(family="SimHei"),#修改背景panel_background=element_blank(),axis_ticks_major_x=element_blank(),axis_ticks_major_y=element_blank(),axis_text=element_text(size=12),axis_title = element_text(size=14),panel_grid_major_x=element_line(color="gray",size=.5),panel_grid_major_y=element_line(color="gray",size=.5),))
Krig_inter_grid.save('test1.jpg',dpi=300)

裁剪后数据可视化

df_grid_geo = gpd.GeoDataFrame(df_grid, geometry=gpd.points_from_xy(df_grid["long"], df_grid["lat"]),crs="EPSG:4326")
#根据之前的js geojson格式的地图文件进行裁剪
js_kde_clip = gpd.clip(df_grid_geo,js)import plotnine
from plotnine import *
plotnine.options.figure_size = (6.4, 4)
Krig_inter_grid = (ggplot() + geom_tile(js_kde_clip,aes(x='long',y='lat',fill='Krig_gaussian'),size=0.1) +geom_map(js,fill='none',color='gray',size=0.6) + scale_fill_cmap(cmap_name='Spectral_r',name='PM2.5',breaks=[30,40,60,80])+scale_x_continuous(breaks=[108,110,112,114,116,118])+scale_y_continuous(breaks=[29,31,33,35])+#labs(title="Map Charts in Python Exercise 02: Map point interpolation",)+#添加文本信息#annotate('text',x=116.5,y=35.3,label="processed map charts with plotnine",ha="left", size=10)+#annotate('text',x=120,y=30.6,label="Visualization by DataCharm",ha="left",size=9)+theme(text=element_text(family="SimHei"),#修改背景panel_background=element_blank(),axis_ticks_major_x=element_blank(),axis_ticks_major_y=element_blank(),axis_text=element_text(size=12),axis_title = element_text(size=14),panel_grid_major_x=element_line(color="gray",size=.5),panel_grid_major_y=element_line(color="gray",size=.5),))
Krig_inter_grid.save('test2.jpg',dpi=300)

参考文献：
https://www.bilibili.com/video/BV1jy4y1m7wU

Python进行克里金插值可视化相关推荐

python：克里金插值
python:克里金插值最近写代码遇到了使用样本数据做克里金插值的事情.于是将Excel保存的[x坐标,y坐标,样本值]数据结合tif数据做了克里金插值,并将代码记录下来. 克里金插值结果: 输入数 ...
【Matlab 克里金】克里金插值
前言早在去年就开始要学 matlab,直到昨天,有明确的需要,和可视化的效果,才定下心来学: 从克里金插值上手,把该会的都捋了一遍. 学习路径克里金插值,有代码,有dacefit函数包,但是第一次 ...
gstat | 空间插值（四）——克里金插值之协同克里金和交叉验证
本篇是介绍克里金插值的第三篇推文,也是最后一篇. 因为前面两篇使用的数据中已知点的样本太少,本篇使用gstat工具包说明文档中的数据集.该数据集来自sp工具包. 实际上,gstat工具包的方法对sp和 ...
gstat | 空间插值（三）——克里金插值之泛克里金和简单克里金
本篇接着上篇继续介绍克里金插值.首先加载相关工具包和上篇使用的示例数据: library(gstat) library(sf) library(tidyverse) library(readxl) l ...
gstat | 空间插值（二）——克里金插值之普通克里金
说明:昨天的推文误把可吸入颗粒物当作PM2.5,实应该为PM10,这里修正后重发. 从本篇开始计划分三篇介绍克里金插值.与反距离权重插值不同,克里金插值是无偏估计,其中也涉及到模型估计.本篇先对普通克 ...
GNSS速度场简易MATLAB克里金插值
引言由于GPS观测点分布离散且不均匀,在进行应变计算和分析前,需要对速度场进行插值,获得均匀分布的速度场.一般采用Kriging 法估计在均匀网格节点上的速度值,需要下载克里金MATLAB工具箱 ...
基于主动学习和克里金插值的空气质量推测
基于主动学习和克里金插值的空气质量推测常慧娟, 於志文, 於志勇, 安琦, 郭斌西北工业大学计算机学院,陕西西安 710072 福州大学数学与计算机科学学院,福建福州 350108 摘要 ...
ArcGIS中使用协同克里金插值（co-kriging interplotation ）对气象数据插值
ArcGIS中如何使用协同克里金插值(co-kriging interplotation )对气象数据插值 ANUSPLIN气象站点数据插值局限性百度搜索ArcGIS 克里金插值搭建梯子搜索Arc ...
克里金插值---MATLAB程序
最近在研究克里金插值,拜读了@lanainluv的笔记,备受启发, 在这里做一些补充,并分享自己的代码,希望对各位有所帮助,有误的地方请批评指正!(有帮助的话点个赞吧~) @lanainluv的原文链 ...

Python进行克里金插值可视化

Python进行克里金插值可视化相关推荐

最新文章

热门文章