Python | 绘制动态放大图与论文散点图

文摘   2024-08-21 00:02   天津  

LXX

读完需要

2
分钟

速读仅需 1 分钟

前言:

最近在上网查资料(摸鱼)的时候看到一幅很有意思的图,本文简单地复现一下。

1

   

复现代码:

import numpy as npimport matplotlib.pyplot as pltfrom matplotlib.transforms import Bbox, TransformedBbox, blended_transform_factoryfrom mpl_toolkits.axes_grid1.inset_locator import BboxPatch, BboxConnectorimport imageioimport rasterio
# 从 GeoTIFF 中读取 NDVI 数据tif_file_path = "C:/Users/LXX/Downloads/NDVI.tif" # 请将路径替换为你的 NDVI 影像文件路径
with rasterio.open(tif_file_path) as src: ndvi_data = src.read(1) # 假设 NDVI 数据在第一个波段 ndvi_transform = src.transform
# 获取 NDVI 图像的范围height, width = ndvi_data.shapedistance = np.arange(width)depth = np.arange(height)
# 创建缩放效果的函数def zoom_effect(ax1, ax2, xmin, xmax, **kwargs): trans1 = blended_transform_factory(ax1.transData, ax1.transAxes) trans2 = blended_transform_factory(ax2.transData, ax2.transAxes) bbox = Bbox.from_extents(xmin, 0, xmax, 1) mybbox1 = TransformedBbox(bbox, trans1) mybbox2 = TransformedBbox(bbox, trans2)
c1 = BboxConnector(mybbox1, mybbox2, loc1=3, loc2=3, **kwargs) c1.set_clip_on(False) c2 = BboxConnector(mybbox1, mybbox2, loc1=4, loc2=4, **kwargs) c2.set_clip_on(False)
bbox_patch1 = BboxPatch(mybbox1, facecolor='none', edgecolor='red', linewidth=1.5) bbox_patch2 = BboxPatch(mybbox2, facecolor='none', edgecolor='red', linewidth=1.5)
ax1.add_patch(bbox_patch1) ax2.add_patch(bbox_patch2) ax1.add_patch(c1) ax1.add_patch(c2)
# 创建主轴和缩放区域的函数def create_zoomed_plot(frame_num): xmin = frame_num * 10 xmax = min((frame_num + 10) * 10 + 150, width)
fig, (ax1, ax2) = plt.subplots(nrows=2, figsize=(8, 10), gridspec_kw={'height_ratios': [5, 1]}) mesh1 = ax1.imshow(ndvi_data, aspect='auto', cmap='jet', origin='lower', extent=[0, width, 0, height]) ax1.set_ylabel('Pixel Index (height)') ax1.set_title('NDVI Zoomed Analysis') cbar1 = fig.colorbar(mesh1, ax=ax1, orientation='vertical', pad=0.01) cbar1.set_label('NDVI')
mesh2 = ax2.imshow(ndvi_data[:, xmin:xmax], aspect='auto', cmap='jet', origin='lower', extent=[xmin, xmax, 0, height]) ax2.set_ylabel('Pixel Index (height)') ax2.set_xlabel('Pixel Index (width)')
zoom_effect(ax1, ax2, xmin, xmax, edgecolor='red', linewidth=1.5)
fig.tight_layout() plt.savefig(f"frame_{frame_num:02d}.png", bbox_inches='tight') plt.close(fig)
# 保存多个帧for i in range(10): create_zoomed_plot(i)
# 使用 imageio 生成 GIFwith imageio.get_writer('ndvi_movie_zoomed.gif', mode='I', fps=1) as writer: for i in range(10): filename = f"frame_{i:02d}.png" image = imageio.imread(filename) writer.append_data(image)
print("GIF 生成完成")


3

   

绘制散点图

散点图(Scatter plot)是一种统计图表,一种数据可视化工具,用于显示两个变量之间的关系。 它通过在坐标平面上绘制数据点来展示变量之间的关联程度。


散点图可以用来帮助我们观察两个变量之间的关系,可以帮助我们看到一些隐藏的问题:
变量之间的趋势:散点图可以告诉我们两个变量之间是如何变化的。如果数据点向右上方倾斜,那么意味着一个变量增加时,另一个变量也会增加。如果数据点向右下方倾斜,那么意味着一个变量增加时,另一个变量会减少。
相关性:散点图可以显示出两个变量之间的相关性。如果数据点在图上聚集成一条直线或者有明显的形状,那么这两个变量可能存在较强的相关关系。如果数据点分散在图上,没有明显的形状或趋势,那么这两个变量可能没有明显的相关关系。
异常值的检测:散点图可以帮助我们发现数据中的异常值。异常值是与其他数据点明显不同的值,可能是由于测量错误或特殊情况导致的。通过观察散点图,我们可以看到是否有一些离其他点很远的点,这些点可能是异常值。

绘制代码:
import numpy as npimport matplotlib.pyplot as pltfrom scipy import stats
# 创建随机数n = 1000x = np.random.uniform(20, 160, n)y = 0.82 * x + 7.68 + np.random.normal(0, 10, n)
# 计算线性拟合slope, intercept, r_value, p_value, std_err = stats.linregress(x, y)fit_line = slope * x + intercept
# 计算均方误差 (MSE)、根均方误差 (RMSE) 和平均绝对误差 (MAE)mse = np.mean((y - fit_line) ** 2)rmse = np.sqrt(mse)mae = np.mean(np.abs(y - fit_line))
# 绘制带拟合线的散点图plt.figure(figsize=(8, 6))scatter = plt.scatter(x, y,marker='.',c=y, cmap='magma', alpha=1, label='采样点')plt.plot(x, fit_line, color='red', label=f'y={slope:.2f}x+{intercept:.2f}')plt.plot([20, 180], [20, 180], '--', color='grey', label='y=x')
plt.rcParams['font.sans-serif'] = ['STSong']plt.rcParams['axes.unicode_minus'] = False# 计算95%置信区间#confidence_interval = 1.96 * std_err#upper_fit_line = fit_line + confidence_interval#lower_fit_line = fit_line - confidence_interval#plt.plot(x, upper_fit_line, 'k-', alpha=0.6)#plt.plot(x, lower_fit_line, 'k-', alpha=0.6)
# 添加注释plt.xlabel('x')plt.ylabel('y')plt.legend(loc='upper left')plt.text(25, 125, f'$R^2$={r_value**2:.2f}, $P$<0.001\nMAE={mae:.2f}\nRMSE={rmse:.2f}\nMAPE={np.mean(np.abs((y - fit_line) / y)):.2f}', fontsize=12)
plt.xlim(20, 180)plt.ylim(20, 180)plt.show()


遥感小屋
分享遥感相关文章、代码,大家一起交流,互帮互助
 最新文章