海洋数据处理软件 -海洋数据统计

$ O8 q# e$ E/ }5 h! k# j

气象海洋领域，常常会涉及到大规模的数据处理，比如高时空分辨率的模式数据、雷达数据和卫星数据等。

$ L% N7 K) e7 ^; G/ X7 q

9 h& b- m) l8 \6 ^3 t) s$ j [( o

大部分情况，我们一般只会涉及小规模数据处理，对计算效率并不会太过追求。但是当数据量变大时，低效的数据处理所耗费的时间非常明显，因此高效的数据处理方式尤为重要。

( F8 E& ?7 V3 E5 x4 e

本篇以拟合一个高维数据的正太分布参数为例，介绍如何使用xarray+dask加速数据处理。

) t6 f% j s- I. @* g

数据维度X[time, lev, lat, lon]，需要对三维空间每一点，沿着时间维度做正太分布拟合(正太分布拟合只是作为例子，这里可以定义你需要的操作函数)。

其中几个关键点解释一下：

) k. [: ]8 \3 V

（1）首先定义拟合正太分布的函数

2 e G) r3 q' l# w4 ?) Y! Z; a

def norm_fit(data): " H6 T2 f& b. y; O loc, scale = st.norm.fit(data) return np.array([loc, scale])

这里需要注意的是，拟合的函数，输出参数，需要打包为一个数组。并且它的维度需要和后面aplly_ufun定义的输出维度一致。

& A9 P3 x9 w% U3 k1 Z

（2）xarray分块

, B, \4 X8 m8 j" L: }) y2 W. `

x = xr.DataArray(rs.randn(500, 20, 1500, 1500), dims=["time", "lev", "lat","lon"]) x = x.chunk({"lev":5, "lat":100,"lon":100})

( d* W; ?% |+ w9 b, L' V" k

xarray结合dask可以将一个大型数组分成一个个数据块（chunk），需要注意的是我们需要沿着时间维操作，拟合需要整个时间维度的数据，因此时间维time不能分块，只能对其他维度分块。

（3）xarray.apply_ufunc函数

result = xr.apply_ufunc(norm_fit, x, " K0 \9 _- {" z7 F/ K2 } input_core_dims=[["time"]], " x+ ^; \- o c$ ^0 |! r! t2 x7 W output_core_dims=[["paras"]], / C3 p* W2 }. i/ X# L dask="parallelized", * i( y6 F# m2 s3 @# ?$ [ output_dtypes=[np.float], 3 q! k# q$ O3 @7 H1 d* Q$ C dask_gufunc_kwargs={output_sizes:{"paras":2}} )

apply_ufunc函数具体可以参考官网教程，这里只说使用时遇到的困难，即如何定义输出维度：输出维度是用output_core_dim定义的，将输出的拟合参数（期望和标准差）定义为paras维度，维度的大小为2，通过output_sizes参数设置。这样我们输入[time, lev, lat, lon]的数据，在每个空间点对时间维度拟合之后，输出的数据为[lev, lat, lon, paras]。（PS: 这里感谢深雨露大佬的指点）

4 ]$ x u# r! `+ j' F% o

以下是全部代码：

from scipy import stats as st / Y; e0 b: D, v# G8 f import xarray as xr 7 W- N5 U/ x1 n9 @ b import dask # V! V g0 b. m! U& Z import numpy as np / V+ f, x" U( {7 H0 d from dask.diagnostics import ProgressBar " o: L9 T- a3 [& A$ p: w def norm_fit(data): loc, scale = st.norm.fit(data) return np.array([loc, scale]) rs = np.random.RandomState(0) " [! f: a, s, [6 l' ^- ` x = xr.DataArray(rs.randn(500, 20, 1500, 1500), dims=["time", "lev", "lat","lon"]) & D1 b0 c, @0 j' r; C3 X2 s4 W x = x.chunk({"lev":5, "lat":100,"lon":100}) ) Z A9 P( D% ^# x+ h* p" X #使用apply_ufunc计算，并用dask的并行计算 result = xr.apply_ufunc(norm_fit, x, input_core_dims=[["time"]], output_core_dims=[["paras"]], # R- N8 D+ q6 \9 ^0 c6 h dask="parallelized", ) _8 j5 g! p! a3 D6 v& C/ f- r output_dtypes=[np.float], ) c( A" w5 r$ \+ A6 n dask_gufunc_kwargs={output_sizes:{"paras":2}} ) 6 M- w6 Q. O/ A0 G8 ^ #compute进行真实的计算并显示进度 0 `& r& S K( q: P% K# Q; f with ProgressBar(): result = results.compute() ! r2 d3 }0 O) Y/ h+ C7 u #结果冲命名保存到nc文件 result = result.rename("norm_paras") & O6 a) Q( Z2 A0 V result = result.to_netcdf("norm_fit_paras.nc",compute=False) with ProgressBar(): % I& A4 e1 o9 `3 H* G' {' k result.compute()

+ |# d/ K) k1 |! j, D1 t

转自：气海同途

关注【Ai尚研修科研技术平台】公众号，查看更多课程安排及免费数据资料

【大气，生态、遥感、水文水资源、地下水环境、统计语言、生信等】多领域视频教程​mp.weixin.qq.com/s?__biz=MzUyNzczMTI4Mg==&mid=2247576733&idx=3&sn=8f147b34d9652c454b954dcc2952ae38&chksm=fa78f3a0cd0f7ab630d262dbe674cd24ee8c06b2b0c12ab029b18a702e692dd6d8c48e0783e6&token=1162449722&lang=zh_CN&scene=21#wechat_redirect

0 }( j/ B2 |* x2 I% y

推荐：

1、Python语言在地球科学领域中的应用实践技术应用

2、【夏令营】针对课题组人员AI培养计划:“开启AI科研之路”

3、Python在气象与海洋中的实践技术应用精品课程

1 ^: `3 E5 X" T9 e$ l6 U

4、Python在WRF模型自动化运行及前后处理中的实践技术应用

5、全套Python机器学习核心技术与案例分析实践应用视频

# D4 z Q. e/ k I

6、全套区域高精度地学模拟-WRF气象建模、多案例应用与精美制图精品课程

! d9 P9 q1 i; I% _ v

7、WRF DA资料同化系统理论、运行与与变分、混合同化新方法技术应用视频教程

! M9 N* G& w8 s/ U' Y! r7 d. d # ^* T' ^( k% t* w8 C: K7 y ) O$ @, U6 C: K0 `