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

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

1 A5 J! ^! \) e) Y7 Q6 h

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

! ~' f" d7 F E/ k. [0 `0 A

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

( v/ ^# L# x$ K: ~# X0 Z, P- d( j * G9 P2 l) w8 q0 D% z+ e& F# F# m

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

4 K& M$ o# B8 M; \

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

: @" M: H7 M+ R/ j

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

def norm_fit(data): 8 ?* J/ g# C% e& ? loc, scale = st.norm.fit(data) - R4 n. P7 }2 n+ M# X& T* M" \ return np.array([loc, scale])

6 c3 _# o0 U* ^3 f

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

1 }4 Y; M. i( O; M7 U. A8 S B

（2）xarray分块

; n# n* a0 K. U: S- Z- f8 I

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

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

（3）xarray.apply_ufunc函数

result = xr.apply_ufunc(norm_fit, 0 w3 G4 }0 I% G0 Q8 { H1 O x, input_core_dims=[["time"]], output_core_dims=[["paras"]], ; o5 ?. W: [* w$ z dask="parallelized", output_dtypes=[np.float], dask_gufunc_kwargs={output_sizes:{"paras":2}} & w7 e& c5 G$ t: M )

5 F3 r; z, j. Q D9 p) V O0 ^; ^

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

6 }1 u: b% g* [' k0 E. Y

以下是全部代码：

from scipy import stats as st import xarray as xr import dask 9 h1 D& Y$ u2 @2 H4 p import numpy as np ! g+ c1 i' z5 |' r: `; E4 q: U) m from dask.diagnostics import ProgressBar def norm_fit(data): 6 i% a6 J* V2 H7 q/ g$ I8 @ loc, scale = st.norm.fit(data) return np.array([loc, scale]) - H/ y2 w, E! S" F/ U1 L rs = np.random.RandomState(0) 1 c. T7 |) q4 B \& b+ E8 o x = xr.DataArray(rs.randn(500, 20, 1500, 1500), dims=["time", "lev", "lat","lon"]) x = x.chunk({"lev":5, "lat":100,"lon":100}) ' l& T1 U6 o% ]$ A( f/ H& m #使用apply_ufunc计算，并用dask的并行计算 * a0 q4 u3 F" @7 L2 i7 U7 O result = xr.apply_ufunc(norm_fit, x, & ?' O3 E# E( c3 } input_core_dims=[["time"]], $ b7 J4 x4 V2 X/ z1 p5 z4 G output_core_dims=[["paras"]], - x% J' j V) E5 l1 A dask="parallelized", output_dtypes=[np.float], dask_gufunc_kwargs={output_sizes:{"paras":2}} + L+ T8 y- a8 ?5 t$ S9 D7 O ) % G! V. {+ o4 v: `7 h" j #compute进行真实的计算并显示进度 $ y/ T8 E7 Q* i with ProgressBar(): result = results.compute() ' D3 ^' {8 w! x4 S+ F #结果冲命名保存到nc文件 result = result.rename("norm_paras") % I9 C, g* N! @+ E" p. Z result = result.to_netcdf("norm_fit_paras.nc",compute=False) 0 r" ~ \2 ^3 b- [& e) H6 [( m with ProgressBar(): result.compute()

+ u: W, w5 {) f

转自：气海同途

: D* }! E$ Y: p0 U. X% z0 P1 A9 O

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

3 z1 A+ I1 Z* U/ Y) l1 e

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

推荐：

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

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

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

$ L' K& ~ ]$ m5 d

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

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

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

/ k {& p s' N& ~

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

3 M7 C: Q+ P% I# P ! N2 W0 R! s2 l' e, m/ P1 S