数值型检验指标
<p>[TOC]</p>
<pre><code class="language-python">%matplotlib inline
%load_ext autoreload
%autoreload 2
import meteva.method as mem
import numpy as np</code></pre>
<p>多分类预报分为两种情况：</p>
<ol>
<li>预报的对象本身不是数，而是类别，比如降水相态的预报和观测量可能是无降水、降雨、降雪、霰、冻雨等，通常我们可以采用整数对不同的类别进行编码。</li>
<li>预报的对象本身是实数，但人为根据数值大小划分为多个等级，从而构成类别，比如根据24小时降水量，记[0.1,10),[10,25),[25,50),[50,100),[100,250),[250,∞）(单位mm/24小时）为小雨、中雨、大雨、暴雨、大暴雨和特大暴雨。 </li>
</ol>
<p>对于上两种情况的预报检验，本程序尽量提供统一函数实现，以提高函数使用的便利性
多分类预报的列联表可以完整的展示预报的性能，根据它可以计算各类统计指标。为了方便展示不同指标的检验方法的内涵，以24小时降水的分类检验为例，一下展示分级降水预报的列联表：</p>
<table>
<thead>
<tr>
<th style="text-align: left;">样本数 /观测<br> ——————— <br>预报\</th>
<th style="text-align: left;">无降水(g0)</th>
<th style="text-align: left;">小雨（g1）</th>
<th style="text-align: left;">中雨(g2)</th>
<th style="text-align: left;">大雨(g3)</th>
<th style="text-align: left;">暴雨(g4)</th>
<th style="text-align: left;">大暴雨(g5)</th>
<th style="text-align: left;">特大暴雨(g6)</th>
<th style="text-align: left;">合计</th>
</tr>
</thead>
<tbody>
<tr>
<td style="text-align: left;">无降水(g0)</td>
<td style="text-align: left;">N_00</td>
<td style="text-align: left;">N_01</td>
<td style="text-align: left;">N_02</td>
<td style="text-align: left;">N_03</td>
<td style="text-align: left;">N_04</td>
<td style="text-align: left;">N_05</td>
<td style="text-align: left;">N_06</td>
<td style="text-align: left;">N_0S</td>
</tr>
<tr>
<td style="text-align: left;">小雨(g1)</td>
<td style="text-align: left;">N_10</td>
<td style="text-align: left;">N_11</td>
<td style="text-align: left;">N_12</td>
<td style="text-align: left;">N_13</td>
<td style="text-align: left;">N_14</td>
<td style="text-align: left;">N_15</td>
<td style="text-align: left;">N_16</td>
<td style="text-align: left;">N_1S</td>
</tr>
<tr>
<td style="text-align: left;">中雨(g2)</td>
<td style="text-align: left;">N_20</td>
<td style="text-align: left;">N_21</td>
<td style="text-align: left;">N_22</td>
<td style="text-align: left;">N_23</td>
<td style="text-align: left;">N_24</td>
<td style="text-align: left;">N_25</td>
<td style="text-align: left;">N_26</td>
<td style="text-align: left;">N_2S</td>
</tr>
<tr>
<td style="text-align: left;">大雨(g3)</td>
<td style="text-align: left;">N_30</td>
<td style="text-align: left;">N_31</td>
<td style="text-align: left;">N_32</td>
<td style="text-align: left;">N_33</td>
<td style="text-align: left;">N_34</td>
<td style="text-align: left;">N_35</td>
<td style="text-align: left;">N_36</td>
<td style="text-align: left;">N_3S</td>
</tr>
<tr>
<td style="text-align: left;">暴雨(g4)</td>
<td style="text-align: left;">N_40</td>
<td style="text-align: left;">N_41</td>
<td style="text-align: left;">N_42</td>
<td style="text-align: left;">N_43</td>
<td style="text-align: left;">N_44</td>
<td style="text-align: left;">N_45</td>
<td style="text-align: left;">N_46</td>
<td style="text-align: left;">N_4S</td>
</tr>
<tr>
<td style="text-align: left;">大暴雨(g5)</td>
<td style="text-align: left;">N_50</td>
<td style="text-align: left;">N_51</td>
<td style="text-align: left;">N_52</td>
<td style="text-align: left;">N_53</td>
<td style="text-align: left;">N_54</td>
<td style="text-align: left;">N_55</td>
<td style="text-align: left;">N_56</td>
<td style="text-align: left;">N_5S</td>
</tr>
<tr>
<td style="text-align: left;">特大暴雨(g6)</td>
<td style="text-align: left;">N_60</td>
<td style="text-align: left;">N_61</td>
<td style="text-align: left;">N_62</td>
<td style="text-align: left;">N_63</td>
<td style="text-align: left;">N_64</td>
<td style="text-align: left;">N_65</td>
<td style="text-align: left;">N_66</td>
<td style="text-align: left;">N_6S</td>
</tr>
<tr>
<td style="text-align: left;">合计</td>
<td style="text-align: left;">N_S0</td>
<td style="text-align: left;">N_S1</td>
<td style="text-align: left;">N_S2</td>
<td style="text-align: left;">N_S3</td>
<td style="text-align: left;">N_S4</td>
<td style="text-align: left;">N_S5</td>
<td style="text-align: left;">N_S6</td>
<td style="text-align: left;">N_total</td>
</tr>
</tbody>
</table>
<p><strong><em>通过随机函数生成测试数据，用于后续检验函数调用示例</em></strong></p>
<pre><code class="language-python">ob = np.random.rand(2,20) * 10
fo1 = np.random.rand(2,20) * 10
fo2 = np.random.rand(5,2,20) * 10
ob_int = ob.astype(np.int8)
fo_int1 = fo1.astype(np.int8)
fo_int2 = fo2.astype(np.int8)
grade_list = [3,5]</code></pre>
<pre><code class="language-python">ob_int</code></pre>
<pre><code>array([[9, 5, 3, 3, 8, 7, 1, 4, 8, 4, 9, 8, 4, 6, 6, 5, 1, 0, 7, 7],
[1, 9, 8, 0, 1, 1, 2, 4, 9, 3, 3, 7, 2, 7, 0, 6, 4, 6, 5, 2]],
dtype=int8)</code></pre>
<pre><code class="language-python">fo_int1</code></pre>
<pre><code>array([[0, 5, 1, 0, 6, 1, 9, 4, 5, 2, 4, 6, 8, 6, 6, 8, 6, 3, 5, 8],
[8, 6, 0, 7, 9, 9, 1, 1, 9, 2, 7, 7, 1, 7, 1, 4, 3, 0, 7, 0]],
dtype=int8)</code></pre>
<pre><code class="language-python">fo_int2</code></pre>
<pre><code>array([[[5, 6, 0, 3, 5, 2, 1, 1, 9, 9, 6, 8, 3, 5, 9, 3, 1, 5, 9, 4],
[7, 0, 9, 0, 2, 9, 6, 2, 9, 9, 6, 8, 6, 0, 8, 1, 9, 7, 0, 6]],
[[7, 3, 0, 6, 7, 7, 5, 0, 6, 9, 9, 0, 8, 1, 5, 0, 1, 2, 3, 8],
[3, 7, 2, 0, 0, 9, 9, 0, 5, 6, 3, 0, 0, 3, 4, 1, 3, 7, 8, 3]],
[[1, 9, 1, 7, 6, 2, 5, 9, 0, 0, 8, 4, 3, 6, 1, 4, 8, 1, 3, 0],
[8, 5, 6, 5, 3, 5, 4, 1, 2, 8, 9, 9, 3, 1, 2, 2, 0, 2, 1, 8]],
[[9, 1, 0, 4, 9, 2, 9, 1, 2, 2, 5, 3, 8, 3, 9, 1, 8, 1, 0, 1],
[5, 4, 5, 9, 8, 5, 0, 4, 0, 2, 6, 6, 2, 6, 9, 4, 2, 0, 4, 4]],
[[0, 1, 9, 3, 2, 1, 4, 4, 5, 2, 5, 1, 6, 4, 6, 4, 6, 5, 6, 0],
[6, 6, 0, 2, 8, 1, 1, 8, 1, 5, 0, 3, 3, 7, 2, 7, 1, 7, 6, 7]]],
dtype=int8)</code></pre>
<h1>准确率</h1>
<p><strong><font face="黑体" color=blue size = 5> accuracy(ob, fo, grade_list=None)</font></strong><br />
统计观测和预报分类一致的样本数占比。 以上述列联表为例，accuracy = (N_00 + N_11 + N_22 + N_33 + N_44 + N_55 + N_66)/N_total</p>
<table>
<thead>
<tr>
<th style="text-align: left;">参数</th>
<th style="text-align: left;">说明</th>
</tr>
</thead>
<tbody>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>ob</font></strong></td>
<td style="text-align: left;">实况数据, 任意维numpy数组</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>fo</font></strong></td>
<td style="text-align: left;">fo比Ob.shape多一维或者保持一致，fo.shape低维与ob.shape保持一致</td>
</tr>
<tr>
<td style="text-align: left;"><strong>grade_list</strong></td>
<td style="text-align: left;">如果该参数为None，观测或预报值出现过的值都作为分类标记，如果该参数不为None，它必须是一个从小到大排列的实数，以其中列出的数值划分出的多个区间作为分类标签。对于预报和观测值不为整数的情况，grade_list 不能设置为None</td>
</tr>
<tr>
<td style="text-align: left;"><font face="黑体" color=blue size=5>return</font></td>
<td style="text-align: left;">如果Fo和Ob的shape一致（即只有一种预报），返回结果为实数；如果Fo比Ob高出一维，则返回1维数组，长度等于预报成员数。元素是 0-1的实数，0代表无技巧，最优预报为1</td>
</tr>
</tbody>
</table>
<p><strong>调用示例：</strong> </p>
<pre><code class="language-python">mem.accuracy(ob_int,fo_int1)</code></pre>
<pre><code>0.175</code></pre>
<pre><code class="language-python">mem.accuracy(ob_int,fo_int1,grade_list)</code></pre>
<pre><code>0.5</code></pre>
<pre><code class="language-python">mem.accuracy(ob_int,fo_int2)</code></pre>
<pre><code>array([0.15 , 0.125, 0.025, 0.075, 0.125])</code></pre>
<pre><code class="language-python">mem.accuracy(ob_int,fo_int2,grade_list)</code></pre>
<pre><code>array([0.475, 0.45 , 0.25 , 0.3 , 0.375])</code></pre>
<p>以下将介绍分类预报和分级预报的检验算法，它们和二分类预报的检验算法的不同可以用如下示意图加以说明：<br />
<img src="https://www.showdoc.com.cn/server/api/attachment/visitfile/sign/4aca708a4afbaf852bb773fcddc422ac" alt="" />
<img src="https://www.showdoc.com.cn/server/api/attachment/visitfile/sign/a750cdca4327bb4be38507b41ab52699" alt="" />
<img src="https://www.showdoc.com.cn/server/api/attachment/visitfile/sign/79195332cb78b52d480014eb465d9b8c" alt="" /></p>
<h1>分类预报ETS评分</h1>
<p><strong><font face="黑体" color=blue size = 5> ets_multi(ob, fo,grade_list=None) </font></strong><br />
基于原始数据计算针对每个类别的ets评分：(Hit-Hit_random) /(Hit + Misses+ False alarms - Hit_random)，<br />
Hit_random 为随机预报的命中数， Hit_random = (Hit + Misses)*(Hit + Falase alarms)/ total<br />
基于上述列联表，以大雨降水分类ts评分计算为例，其中 </p>
<ul>
<li>Hit = N_33</li>
<li>Misses = (N_S3 - N_33) </li>
<li>False alarms = (N_3S - N_33) </li>
<li>total = N_total。 </li>
</ul>
<h1>分类预报TS评分</h1>
<p><strong><font face="黑体" color=blue size = 5>ts_multi(ob, fo,grade_list=None)</font></strong><br />
基于原始数据计算针对每个类别的ts评分： Hit /(Hit + Misses+ False alarms) ,其中Hit,Misses 和False alarms 含义同上。 </p>
<h1>分类预报BIAS</h1>
<p><strong><font face="黑体" color=blue size = 5>bias_multi(ob, fo,grade_list=None)</font></strong><br />
基于原始数据计算针对每个类别的bias评分：(Hit + False alarms)/(Hit + Misses)，其中Hit,Misses 和False alarms 含义同上。</p>
<h1>分类预报空报率</h1>
<p><strong><font face="黑体" color=blue size = 5>far_multi(ob, fo,grade_list=None)</font></strong><br />
基于原始数据计算针对每个类别的far评分： False alarms/(Hit + False alarms)，其中Hit,Misses 和False alarms 含义同上。</p>
<h1>分类预报漏报率</h1>
<p><strong><font face="黑体" color=blue size = 5>mr_multi(ob, fo,grade_list=None)</font></strong><br />
基于原始数据计算针对每个类别的mr评分： Misses/(Hit + Misses)，其中Hit,Misses 和False alarms 含义同上。</p>
<h1>分类预报成功率</h1>
<p><strong><font face="黑体" color=blue size = 5>sr_multi(ob, fo,grade_list=None)</font></strong><br />
基于原始数据计算针对每个类别的sr评分： Hits/(Hits + False alarms)，其中Hit 和False alarms 含义同上。</p>
<h1>分类预报命中率</h1>
<p><strong><font face="黑体" color=blue size = 5>pod_multi(ob, fo,grade_list=None)</font></strong><br />
基于原始数据计算针对每个类别的pod评分： Hits/(Hits + Misses)，其中Hit 和Misses 含义同上。</p>
<h1>分类预报报空率</h1>
<p><strong><font face="黑体" color=blue size = 5>pofd_multi(ob, fo,grade_list=None)</font></strong><br />
基于原始数据计算针对每个类别的pod评分： False alarms/(False alarms + Correct negatives)，其中False alarms含义同上
Correct negatives = total - Hit -Misses - False alarms。</p>
<table>
<thead>
<tr>
<th style="text-align: left;">参数</th>
<th style="text-align: left;">说明</th>
</tr>
</thead>
<tbody>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>ob</font></strong></td>
<td style="text-align: left;">实况数据, 任意维numpy数组</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>fo</font></strong></td>
<td style="text-align: left;">fo比Ob.shape多一维或者保持一致，fo.shape低维与ob.shape保持一致</td>
</tr>
<tr>
<td style="text-align: left;"><strong>grade_list</strong></td>
<td style="text-align: left;">如果该参数为None，观测或预报值出现过的值都作为分类标记，并由此确定类别数。如果该参数不为None，它必须是一个从小到大排列的实数，以其中列出的数值划分出的多个区间作为分类标签，例如当grade_list = [3,5]，样本将被分为（负无穷,3),[3,5),[5,正无穷）三种类别。对于预报和观测值不为整数的情况，grade_list 不能设置为None</td>
</tr>
<tr>
<td style="text-align: left;"><font face="黑体" color=blue size=5>return</font></td>
<td style="text-align: left;">当fo包含多个成员预报时，返回一个多维numpy数组,shaop= (预报成员数，类别数),否则返回一维数组，长度等于类别数。数组中每个元素为一个模式一个类别的ts评分 。每个元素值为0到1的实数，完美预报对应值为1。</td>
</tr>
</tbody>
</table>
<p> </p>
<p> </p>
<p><em>------------------------------</em></p>
<h1>分级预报ETS评分</h1>
<p><strong><font face="黑体" color=blue size = 5>ets_grade(ob, fo,grade_list=None)</font></strong><br />
《全国智能网格陆面基本要素预报检验办法》中表3.1对应的分级ets评分：(Hit-Hit_random) /(Hit + Misses+ False alarms - Hit_random)，<br />
Hit_random 为随机预报的命中数， Hit_random = (Hit + Misses)*(Hit + Falase alarms)/ total 。<br />
基于上述列联表，以大雨降水分类ts评分计算为例，其中 </p>
<ul>
<li>Hit = N_33</li>
<li>Misses = (N_03 + N_13 + N_23) </li>
<li>False alarms = (N_30 + N_31 + N_32) </li>
<li>total = sum(N_ij),(i = 0,1,2,3, j = 0,1,2,3)。 </li>
</ul>
<h1>分级预报TS评分</h1>
<p><strong><font face="黑体" color=blue size = 5>ts_grade(ob, fo,grade_list=None)</font></strong><br />
《全国智能网格陆面基本要素预报检验办法》中表3.1对应的分级ts评分：Hit /(Hit + Misses+ False alarms)，其中Hit,Misses 和False alarms 含义同上。 </p>
<h1>分级预报BIAS评分</h1>
<p><strong><font face="黑体" color=blue size = 5>bias_grade(ob, fo,grade_list=None)</font></strong><br />
《全国智能网格陆面基本要素预报检验办法》中表3.1对应的分级BIAS评分：(Hit + False alarms)/(Hit + Misses)，其中Hit,Misses 和False alarms 含义同上。 </p>
<h1>分级预报空报率</h1>
<p><strong><font face="黑体" color=blue size = 5>far_grade(ob, fo,grade_list=None)</font></strong><br />
《全国智能网格陆面基本要素预报检验办法》中表3.1对应的分级空报率：False alarms/(Hit + False alarms)，其中Hit,Misses 和False alarms 含义同上。 </p>
<h1>分级预报漏报率</h1>
<p><strong><font face="黑体" color=blue size = 5>mr_grade(ob, fo,grade_list=None)</font></strong><br />
《全国智能网格陆面基本要素预报检验办法》中表3.1对应的分级漏报率：Misses/(Hit + Misses)，其中Hit,Misses 和False alarms 含义同上。 </p>
<h1>分级预报成功率</h1>
<p><strong><font face="黑体" color=blue size = 5>sr_grade(ob, fo,grade_list=None)</font></strong><br />
基于原始数据计算针对每个级别的sr评分： Hits/(Hits + False alarms)，其中Hit 和False alarms 含义同上。</p>
<h1>分级预报命中率</h1>
<p><strong><font face="黑体" color=blue size = 5>pod_grade(ob, fo,grade_list=None)</font></strong><br />
基于原始数据计算针对每个级别的pod评分： Hits/(Hits + Misses)，其中Hit 和Misses 含义同上。</p>
<h1>分级预报报空率</h1>
<p><strong><font face="黑体" color=blue size = 5>pofd_grade(ob, fo,grade_list=None)</font></strong><br />
基于原始数据计算针对每个级别的pod评分： False alarms/(False alarms + Correct negatives)，其中False alarms含义同上
Correct negatives = total - Hit -Misses - False alarms。</p>
<table>
<thead>
<tr>
<th style="text-align: left;">参数</th>
<th style="text-align: left;">说明</th>
</tr>
</thead>
<tbody>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>ob</font></strong></td>
<td style="text-align: left;">实况数据, 任意维numpy数组</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>fo</font></strong></td>
<td style="text-align: left;">fo比Ob.shape多一维或者保持一致，fo.shape低维与ob.shape保持一致</td>
</tr>
<tr>
<td style="text-align: left;"><strong>grade_list</strong></td>
<td style="text-align: left;">必须是一个从小到大排列的实数，以其中列出的数值划分出的多个区间作为分类标签，例如当grade_list = [3,5]，样本将被分为（负无穷,3),[3,5),[5,正无穷）三种类别</td>
</tr>
<tr>
<td style="text-align: left;"><font face="黑体" color=blue size=5>return</font></td>
<td style="text-align: left;">当fo包含多个成员预报时，返回一个多维numpy数组,shaop= (预报成员数，类别数),否则返回一维数组，长度等于类别数。数组中每个元素为一个模式一个类别的ts评分 。每个元素值为0到1的实数，完美预报对应值为1。</td>
</tr>
</tbody>
</table>
<p>对于《全国智能网格陆面基本要素预报检验办法》中表3.2对应的累加降水量级TS评分，应该采用本程序提供的mem.ts方法来计算，实际上就是将降水以某个阈值分割为两种情况，一种是达到阈值记为事件发生，另一种就是为达到阈值记为事件未发生。基于上述列联表，以大雨降水累加降水量级TS评分为例，ts = Hit /(Hit + Misses+ False alarms),其中 hit = sum(N_ij)，（i = 3,4,5,6, j = 3,4,5,6）；Misses = sum(N_ij), （i = 0,1,2,j = 3,4,5,6）；false alarms = sum(N_ij),（ i= 3,4,5,6,j = 0,1,2）。</p>
<p><strong>调用示例：</strong> </p>
<pre><code class="language-python">#grade_list = None,采用观测预报中出现过的0,1,2,3,4,5,6,7,8,9分别被作为类别标签，
mem.ts_multi(ob_int,fo_int1) </code></pre>
<pre><code>array([0. , 0. , 0. , 0. , 0.14285714,
0.2 , 0.25 , 0.25 , 0. , 0.14285714])</code></pre>
<pre><code class="language-python">#grade_list = [3,5],将预报和观测分为（负无穷,3),[3,5),[5,正无穷）三种类别。
mem.ts_multi(ob_int,fo_int1,grade_list)</code></pre>
<pre><code>array([0.2 , 0.16666667, 0.5 ])</code></pre>
<pre><code class="language-python">#grade_list = [3,5],将预报和观测分为（负无穷,3),[3,5),[5,正无穷）三种类别。
mem.ts_multi(ob_int,fo_int2,grade_list)</code></pre>
<pre><code>array([[0.21052632, 0.18181818, 0.41935484],
[0.25 , 0.13333333, 0.40740741],
[0.08 , 0.06666667, 0.23333333],
[0.125 , 0.13333333, 0.24137931],
[0.18181818, 0.14285714, 0.31034483]])</code></pre>
<pre><code class="language-python">#分类ts评分算法的底层仍然是二分类预报的Ts评分算法，只是正负样本的判断标准的差异。
#在本程序库的二分类预报的ts评分函数中也可以接受grade_list参数，通过下面的例子可以了解其中差异
#在ts函数中，grade_list =[3,5]对应计算了两个评分值
#第一个评分代表将预报观测值分为（负无穷，3），[3,正无穷）
#第二个评分代表将预报观测值分为（负无穷，5），[5，正无穷）
#虽然grade_list和上一个例子中是一样的，但返回结果的size却少了1
mem.ts(ob_int,fo_int2,grade_list)</code></pre>
<pre><code>array([[0.58333333, 0.41935484],
[0.57142857, 0.40740741],
[0.39473684, 0.23333333],
[0.43243243, 0.24137931],
[0.5 , 0.31034483]])</code></pre>
<pre><code class="language-python">mem.ets_multi(ob_int,fo_int1,grade_list) #分类预报ETS评分</code></pre>
<pre><code>array([0.02587519, 0.08045977, 0.17647059])</code></pre>
<pre><code class="language-python">mem.bias_multi(ob_int,fo_int1,grade_list) #分类预报bias</code></pre>
<pre><code>array([1.18181818, 0.55555556, 1.1 ])</code></pre>
<pre><code class="language-python">mem.far_multi(ob_int,fo_int1,grade_list) #分类预报空报率</code></pre>
<pre><code>array([0.69230769, 0.6 , 0.36363636])</code></pre>
<pre><code class="language-python">mem.mr_multi(ob_int,fo_int1,grade_list) #分类预报漏报率</code></pre>
<pre><code>array([0.63636364, 0.77777778, 0.3 ])</code></pre>
<pre><code class="language-python">mem.ts_grade(ob_int,fo_int1,grade_list) #分级预报TS评分</code></pre>
<pre><code>array([1. , 0.25, 0.5 ])</code></pre>
<pre><code class="language-python">mem.ets_grade(ob_int,fo_int1,grade_list) #分级预报ETS评分</code></pre>
<pre><code>array([9.99999000e+05, 4.00000000e-02, 1.76470588e-01])</code></pre>
<pre><code class="language-python">mem.bias_grade(ob_int,fo_int1,grade_list) #分级预报bias评分</code></pre>
<pre><code>array([1. , 0.42857143, 1.1 ])</code></pre>
<pre><code class="language-python">mem.far_grade(ob_int,fo_int1,grade_list) #分级预报空报率</code></pre>
<pre><code>array([0. , 0.33333333, 0.36363636])</code></pre>
<pre><code class="language-python">mem.mr_grade(ob_int,fo_int1,grade_list) #分级预报漏报率</code></pre>
<pre><code>array([0. , 0.71428571, 0.3 ])</code></pre>
<h1>hss评分</h1>
<p><strong><font face="黑体" color=blue size = 5>hss(ob,fo,grade_list = None)</font></strong> </p>
<p>Heidke skill score，统计准确率相对于随机预报的技巧</p>
<table>
<thead>
<tr>
<th style="text-align: left;">参数</th>
<th style="text-align: left;">说明</th>
</tr>
</thead>
<tbody>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>ob</font></strong></td>
<td style="text-align: left;">实况数据, 任意维numpy数组</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>fo</font></strong></td>
<td style="text-align: left;">fo比Ob.shape多一维或者保持一致，fo.shape低维与ob.shape保持一致</td>
</tr>
<tr>
<td style="text-align: left;"><strong>grade_list</strong></td>
<td style="text-align: left;">如果该参数为None，观测或预报值出现过的值都作为分类标记，如果该参数不为None，它必须是一个从小到大排列的实数，以其中列出的数值划分出的多个区间作为分类标签。对于预报和观测值不为整数的情况，grade_list 不能设置为None</td>
</tr>
<tr>
<td style="text-align: left;"><font face="黑体" color=blue size=5>return</font></td>
<td style="text-align: left;">如果Fo和Ob的shape一致（即只有一种预报），返回结果为实数；如果Fo比Ob高出一维，则返回1维数组，长度等于预报成员数。-1到1的实数，0代表无技巧，最优预报为1</td>
</tr>
</tbody>
</table>
<p><strong>调用示例：</strong> </p>
<pre><code class="language-python">mem.hss(ob_int,fo_int1)</code></pre>
<pre><code>0.0807799442896936</code></pre>
<pre><code class="language-python">mem.hss(ob_int,fo_int2)</code></pre>
<pre><code>array([ 0.05882353, 0.03114187, -0.08183079, -0.02493075, 0.02574809])</code></pre>
<pre><code class="language-python">mem.hss(ob_int,fo_int1,grade_list)</code></pre>
<pre><code>0.17695473251028807</code></pre>
<pre><code class="language-python">mem.hss(ob_int,fo_int2,grade_list)</code></pre>
<pre><code>array([ 0.11764706, 0.1321499 , -0.17531832, -0.08527132, 0.01185771])</code></pre>
<h1>hk评分</h1>
<p><strong><font face="黑体" color=blue size = 5>hk(ob,fo,grade_list = None)</font></strong><br />
Hanssen and Kuipers discriminant，统计准确率相对于随机预报的技巧</p>
<table>
<thead>
<tr>
<th style="text-align: left;">参数</th>
<th style="text-align: left;">说明</th>
</tr>
</thead>
<tbody>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>ob</font></strong></td>
<td style="text-align: left;">实况数据, 任意维numpy数组</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>fo</font></strong></td>
<td style="text-align: left;">fo比Ob.shape多一维或者保持一致，fo.shape低维与ob.shape保持一致</td>
</tr>
<tr>
<td style="text-align: left;"><strong>grade_list</strong></td>
<td style="text-align: left;">如果该参数为None，观测或预报值出现过的值都作为分类标记，如果该参数不为None，它必须是一个从小到大排列的实数，以其中列出的数值划分出的多个区间作为分类标签。对于预报和观测值不为整数的情况，grade_list 不能设置为None</td>
</tr>
<tr>
<td style="text-align: left;"><font face="黑体" color=blue size=5>return</font></td>
<td style="text-align: left;">如果Fo和Ob的shape一致（即只有一种预报），返回结果为实数；如果Fo比Ob高出一维，则返回1维数组，长度等于预报成员数。-1到1的实数，0代表无技巧，最优预报为1</td>
</tr>
</tbody>
</table>
<p><strong>调用示例：</strong> </p>
<pre><code class="language-python">mem.hk(ob_int,fo_int1)</code></pre>
<pre><code>0.0808926080892608</code></pre>
<pre><code class="language-python">mem.hk(ob_int,fo_int1,grade_list = [1.5,3.5])</code></pre>
<pre><code>0.06844547563805102</code></pre>
<pre><code class="language-python">mem.hss_yesorno(ob_int,fo_int1,grade_list = [1.5,3.5])</code></pre>
<pre><code>array([-0.16438356, 0.25333333])</code></pre>
<h1>seeps误差</h1>
<p><strong><font face="黑体" color=blue size = 5>seeps(ob, fo, p1=None, threshold=None)</font></strong><br />
SEEPS检验指标，具体方法请参考[1]陈法敬, 陈静. "SEEPS"降水预报检验评分方法在我国降水预报中的应用试验[J]. 气象科技进展, 2015(5).</p>
<table>
<thead>
<tr>
<th style="text-align: left;">参数</th>
<th style="text-align: left;">说明</th>
</tr>
</thead>
<tbody>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>ob</font></strong></td>
<td style="text-align: left;">实况数据, 任意维numpy数组</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>fo</font></strong></td>
<td style="text-align: left;">fo比Ob.shape多一维或者保持一致，fo.shape低维与ob.shape保持一致</td>
</tr>
<tr>
<td style="text-align: left;"><strong>p1</strong></td>
<td style="text-align: left;">无雨（观测降水≤0.2mm）的气候概率，缺省时，程序根据ob计算观测降水≤0.2mm的比例作为该气候概率的代替值</td>
</tr>
<tr>
<td style="text-align: left;"><strong>threshold</strong></td>
<td style="text-align: left;">大雨阈值，此处大雨并非指业务标准中的25mm，而是根据降水气候概率分布统计出的较为显著的降水，缺省时会根据输入数据ob来统计该阈值</td>
</tr>
<tr>
<td style="text-align: left;"><font face="黑体" color=blue size=5>return</font></td>
<td style="text-align: left;">如果Fo和Ob的shape一致（即只有一种预报），返回结果为实数；如果Fo比Ob高出一维，则返回1维数组，长度等于预报成员数。常数预报seeps = 1，最优预报seeps=0</td>
</tr>
</tbody>
</table>
<p><strong>调用示例：</strong> </p>
<pre><code class="language-python">mem.seeps(ob,fo1)</code></pre>
<pre><code>1.0441595441595442</code></pre>
<h1>seeps技巧</h1>
<p><strong><font face="黑体" color=blue size = 5>seeps_skill(ob, fo, p1=None, threshold=None)</font></strong><br />
即 1- seeps</p>
<table>
<thead>
<tr>
<th style="text-align: left;">参数</th>
<th style="text-align: left;">说明</th>
</tr>
</thead>
<tbody>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>ob</font></strong></td>
<td style="text-align: left;">实况数据, 任意维numpy数组</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>fo</font></strong></td>
<td style="text-align: left;">fo比Ob.shape多一维或者保持一致，fo.shape低维与ob.shape保持一致</td>
</tr>
<tr>
<td style="text-align: left;"><strong>p1</strong></td>
<td style="text-align: left;">无雨（观测降水≤0.2mm）的气候概率，缺省时，程序根据ob计算观测降水≤0.2mm的比例作为该气候概率的代替值</td>
</tr>
<tr>
<td style="text-align: left;"><strong>threshold</strong></td>
<td style="text-align: left;">大雨阈值，此处大雨并非指业务标准中的25mm，而是根据降水气候概率分布统计出的较为显著的降水，缺省时会根据输入数据ob来统计该阈值</td>
</tr>
<tr>
<td style="text-align: left;"><font face="黑体" color=blue size=5>return</font></td>
<td style="text-align: left;">如果Fo和Ob的shape一致（即只有一种预报），返回结果为实数；如果Fo比Ob高出一维，则返回1维数组，长度等于预报成员数。常数预报seeps_skill = 0，最优预报seeps_skill=1</td>
</tr>
</tbody>
</table>
<p><strong>调用示例：</strong> </p>
<pre><code class="language-python">mem.seeps_skill(ob,fo1)</code></pre>
<pre><code>-0.04415954415954415</code></pre>
<p>在以上示例中，观测和预报的数据都直接输入到评分函数中进行计算，然而有些情况下待检验的数据太大不能整体存入一个numpy数组中，或者不方便整体存入一个numpy数组中，此时就不能调用上面的方式调用评分函数。为此本程序库设计了一套可分块计算的检验程序。<br />
其检验步骤如下：<br />
<strong><em>步骤1：根据需要将分块数据逐一输入到中间结果计算函数</em></strong><br />
<strong><em>步骤2：将中间结果进行累加或合并</em></strong><br />
<strong><em>步骤3：根据累加或合并的中间结果计算检验指标</em></strong><br />
通常上述计算中步骤1是最耗费计算资源，为了提高效率步骤1也可以采用<strong>并行</strong>的方式执行。此外，步骤1执行的结果也可<strong>输出到文件</strong>中，在后续的检验可以从中读入部分中间结果执行后续步骤，从而可以实现各种方式的分组检验，大大提高检验计算效率。</p>
<h1>总样本数、正确样本数</h1>
<p><strong><font face="黑体" color=blue size = 5>tc(ob, fo, grade_list=None)</font></strong><br />
用来计算accuracy等检验指标的中间量 </p>
<table>
<thead>
<tr>
<th style="text-align: left;">参数</th>
<th style="text-align: left;">说明</th>
</tr>
</thead>
<tbody>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>ob</font></strong></td>
<td style="text-align: left;">实况数据, 任意维numpy数组</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>fo</font></strong></td>
<td style="text-align: left;">fo比Ob.shape多一维或者保持一致，fo.shape低维与ob.shape保持一致</td>
</tr>
<tr>
<td style="text-align: left;"><strong>grade_list</strong></td>
<td style="text-align: left;">如果该参数为None，观测或预报值出现过的值都作为分类标记，并由此确定类别数。如果该参数不为None，它必须是一个从小到大排列的实数，以其中列出的数值划分出的多个区间作为分类标签，例如当grade_list = [3,5]，样本将被分为（负无穷,3),[3,5),[5,正无穷）三种类别。对于预报和观测值不为整数的情况，grade_list 不能设置为None</td>
</tr>
<tr>
<td style="text-align: left;"><font face="黑体" color=blue size=5>return</font></td>
<td style="text-align: left;">如果Fo和Ob的shape一致（即只有一种预报），返回结果为shape= （类别数，2）的数组；如果Fo比Ob高出一维，则shape = (预报成员数，类别数，2）的数组，最后一维的内容其依次为总样本数、正确分类的样本数</td>
</tr>
</tbody>
</table>
<h1>总样本数、正确样本数、 预报观测频率表</h1>
<p><strong><font face="黑体" color=blue size = 5>tcof(ob,fo,grade_list = None)</font></strong><br />
用来计算accuracy等检验指标的中间量 </p>
<table>
<thead>
<tr>
<th style="text-align: left;">参数</th>
<th style="text-align: left;">说明</th>
</tr>
</thead>
<tbody>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>ob</font></strong></td>
<td style="text-align: left;">实况数据, 任意维numpy数组</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>fo</font></strong></td>
<td style="text-align: left;">fo比Ob.shape多一维或者保持一致，fo.shape低维与ob.shape保持一致</td>
</tr>
<tr>
<td style="text-align: left;"><strong>grade_list</strong></td>
<td style="text-align: left;">如果该参数为None，观测或预报值出现过的值都作为分类标记，并由此确定类别数。如果该参数不为None，它必须是一个从小到大排列的实数，以其中列出的数值划分出的多个区间作为分类标签，例如当grade_list = [3,5]，样本将被分为（负无穷,3),[3,5),[5,正无穷）三种类别。对于预报和观测值不为整数的情况，grade_list 不能设置为None</td>
</tr>
<tr>
<td style="text-align: left;"><font face="黑体" color=blue size=5>return</font></td>
<td style="text-align: left;">如果Fo和Ob的shape一致（即只有一种预报），返回结果为shape= （类别数+1，2）的数组；如果Fo比Ob高出一维，则shape = (预报成员数，类别数+1，2）的数组。记为返回值为re，其中re[:,0,0] 为总样本数，re[:,0,1]为正确样本数据、re[:,1:,0]为观测频率表，re[:,1:,1]为预报频率表</td>
</tr>
</tbody>
</table>
<h1>分类检验的命中、空报、漏报、正确否定</h1>
<p><strong><font face="黑体" color=blue size = 5>hfmc_multi(ob, fo, grade_list=None)</font></strong><br />
用来计算常用二分类预报检验指标的中间统计量 </p>
<table>
<thead>
<tr>
<th style="text-align: left;">参数</th>
<th style="text-align: left;">说明</th>
</tr>
</thead>
<tbody>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>ob</font></strong></td>
<td style="text-align: left;">实况数据, 任意维numpy数组</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>fo</font></strong></td>
<td style="text-align: left;">fo比Ob.shape多一维或者保持一致，fo.shape低维与ob.shape保持一致</td>
</tr>
<tr>
<td style="text-align: left;"><strong>grade_list</strong></td>
<td style="text-align: left;">如果该参数为None，观测或预报值出现过的值都作为分类标记，并由此确定类别数。如果该参数不为None，它必须是一个从小到大排列的实数，以其中列出的数值划分出的多个区间作为分类标签，例如当grade_list = [3,5]，样本将被分为（负无穷,3),[3,5),[5,正无穷）三种类别。对于预报和观测值不为整数的情况，grade_list 不能设置为None</td>
</tr>
<tr>
<td style="text-align: left;"><font face="黑体" color=blue size=5>return</font></td>
<td style="text-align: left;">如果Fo和Ob的shape一致（即只有一种预报），返回结果为shape= （类别数，4）的数组；如果Fo比Ob高出一维，则shape = (预报成员数，类别数，4）的数组，最后一维的内容其依次为命中、空报、漏报、正确否定的样本数</td>
</tr>
</tbody>
</table>
<h1>分级检验的命中、空报、漏报、正确否定</h1>
<p><strong><font face="黑体" color=blue size = 5>hfmc_grade(ob, fo, grade_list=None)</font></strong><br />
用来计算常用二分类预报检验指标的中间统计量 </p>
<table>
<thead>
<tr>
<th style="text-align: left;">参数</th>
<th style="text-align: left;">说明</th>
</tr>
</thead>
<tbody>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>ob</font></strong></td>
<td style="text-align: left;">实况数据, 任意维numpy数组</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>fo</font></strong></td>
<td style="text-align: left;">fo比Ob.shape多一维或者保持一致，fo.shape低维与ob.shape保持一致</td>
</tr>
<tr>
<td style="text-align: left;"><strong>grade_list</strong></td>
<td style="text-align: left;">必须是一个从小到大排列的实数，以其中列出的数值划分出的多个区间作为分类标签，例如当grade_list = [3,5]，样本将被分为（负无穷,3),[3,5),[5,正无穷）三种类别</td>
</tr>
<tr>
<td style="text-align: left;"><font face="黑体" color=blue size=5>return</font></td>
<td style="text-align: left;">如果Fo和Ob的shape一致（即只有一种预报），返回结果为shape= （类别数，4）的数组；如果Fo比Ob高出一维，则shape = (预报成员数，类别数，4）的数组，最后一维的内容其依次为命中、空报、漏报、正确否定的样本数</td>
</tr>
</tbody>
</table>
<p><strong>调用示例：</strong> </p>
<pre><code class="language-python">tc_array2 = np.zeros((5,2))
tcof_array2 = np.zeros((5,len(grade_list)+2,2))
hfmc_array2 = np.zeros((5,len(grade_list)+1,4))
hfmc_grade2 = np.zeros((5,len(grade_list)+1,4))
for i in range(2):
ob1 = ob_int[i,:]
fo1 = fo_int2[:,i,:]
tc_array2 += mem.tc(ob1,fo1,grade_list)
tcof_array2 += mem.tcof(ob1,fo1,grade_list)
hfmc_array2 += mem.hfmc_multi(ob1,fo1,grade_list)
hfmc_grade2 += mem.hfmc_grade(ob1,fo1,grade_list)
</code></pre>
<pre><code class="language-python">tc_array1 = np.zeros((2))
tcof_array1 = np.zeros((len(grade_list)+2,2))
hfmc_array1 = np.zeros((len(grade_list)+1,4))
hfmc_grade1 = np.zeros((len(grade_list)+1,4))
for i in range(2):
ob1 = ob_int[i,:]
fo1 = fo_int1[i,:]
tc_array1 += mem.tc(ob1,fo1,grade_list)
tcof_array1 += mem.tcof(ob1,fo1,grade_list)
hfmc_array1 += mem.hfmc_multi(ob1,fo1,grade_list)
hfmc_grade1 += mem.hfmc_grade(ob1,fo1,grade_list)
</code></pre>
<p><strong>以下为根据合并后的中间统计量计算最终检验指标的函数：</strong> </p>
<h1>准确率（并行1）</h1>
<p><strong><font face="黑体" color=blue size = 5>accuracy_tc(tc_array)</font></strong><br />
统计观测和预报分类一致的样本数占比 </p>
<table>
<thead>
<tr>
<th style="text-align: left;">参数</th>
<th style="text-align: left;">说明</th>
</tr>
</thead>
<tbody>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>tc_array</font></strong></td>
<td style="text-align: left;">包含检验中间结果的多维数组，其中最后一维长度为2，分别记录了总样本数、正确分类的样本数，它通常是tc函数的计算结果，或者计算结果的累加</td>
</tr>
<tr>
<td style="text-align: left;"><font face="黑体" color=blue size=5>return</font></td>
<td style="text-align: left;">整数或数组，它比tc_array少了最后一维。其中每个元素为0到1的实数，0代表无技巧，最优值为1</td>
</tr>
</tbody>
</table>
<p><strong>调用示例：</strong> </p>
<pre><code class="language-python">mem.accuracy_tc(tc_array1)</code></pre>
<pre><code>0.5</code></pre>
<pre><code class="language-python">mem.accuracy_tc(tc_array2)</code></pre>
<pre><code>array([0.475, 0.45 , 0.25 , 0.3 , 0.375])</code></pre>
<h1>准确率（并行2）</h1>
<p><strong><font face="黑体" color=blue size = 5>accuracy_tcof(tcof_array)</font></strong><br />
统计观测和预报分类一致的样本数占比</p>
<table>
<thead>
<tr>
<th style="text-align: left;">参数</th>
<th style="text-align: left;">说明</th>
</tr>
</thead>
<tbody>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>tc_array</font></strong></td>
<td style="text-align: left;">包含检验中间结果的多维数组，它通常是tcof函数的计算结果，或者计算结果的累加</td>
</tr>
<tr>
<td style="text-align: left;"><font face="黑体" color=blue size=5>return</font></td>
<td style="text-align: left;">整数或数组，它比tcof_array少了最后一维。其中每个元素为0到1的实数，0代表无技巧，最优值为1</td>
</tr>
</tbody>
</table>
<p><strong>调用示例：</strong> </p>
<pre><code class="language-python">mem.accuracy_tcof(tcof_array1)</code></pre>
<pre><code>0.5</code></pre>
<pre><code class="language-python">mem.accuracy_tcof(tcof_array2)</code></pre>
<pre><code>array([0.475, 0.45 , 0.25 , 0.3 , 0.375])</code></pre>
<h1>hss评分（并行）</h1>
<p><strong><font face="黑体" color=blue size = 5>hss_tcof(tcof_array)</font></strong><br />
Heidke skill score，统计准确率相对于随机预报的技巧</p>
<table>
<thead>
<tr>
<th style="text-align: left;">参数</th>
<th style="text-align: left;">说明</th>
</tr>
</thead>
<tbody>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>tc_array</font></strong></td>
<td style="text-align: left;">包含检验中间结果的多维数组，它通常是tcof函数的计算结果，或者计算结果的累加</td>
</tr>
<tr>
<td style="text-align: left;"><font face="黑体" color=blue size=5>return</font></td>
<td style="text-align: left;">整数或数组，它比tcof_array少了最后一维。其中每个元素为-1到1的实数，0代表无技巧，最优预报为1</td>
</tr>
</tbody>
</table>
<p><strong>调用示例：</strong> </p>
<pre><code class="language-python">mem.hss_tcof(tcof_array1)</code></pre>
<pre><code>0.17695473251028807</code></pre>
<pre><code class="language-python">mem.hss_tcof(tcof_array2)</code></pre>
<pre><code>array([ 0.11764706, 0.1321499 , -0.17531832, -0.08527132, 0.01185771])</code></pre>
<h1>hk评分（并行）</h1>
<p><strong><font face="黑体" color=blue size = 5>hk_tcof(tcof_array)</font></strong><br />
Hanssen and Kuipers discriminant，统计准确率相对于随机预报的技巧</p>
<table>
<thead>
<tr>
<th style="text-align: left;">参数</th>
<th style="text-align: left;">说明</th>
</tr>
</thead>
<tbody>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>tc_array</font></strong></td>
<td style="text-align: left;">包含检验中间结果的多维数组，它通常是tcof函数的计算结果，或者计算结果的累加</td>
</tr>
<tr>
<td style="text-align: left;"><font face="黑体" color=blue size=5>return</font></td>
<td style="text-align: left;">整数或数组，它比tcof_array少了最后一维。其中每个元素为-1到1的实数，0代表无技巧，最优预报为1</td>
</tr>
</tbody>
</table>
<p><strong>调用示例：</strong> </p>
<pre><code class="language-python">mem.hk_tcof(tcof_array1)</code></pre>
<pre><code>0.1723446893787575</code></pre>
<pre><code class="language-python">mem.hk_tcof(tcof_array2)</code></pre>
<pre><code>array([ 0.11222445, 0.13426854, -0.17935872, -0.08817635, 0.01202405])</code></pre>
<p>基于hfmc计算分类预报的ts、ets、bias、far、mr等评分时，其函数调用方法和二分类模块中的方法是一样的。</p>
<pre><code class="language-python">mem.ts_hfmc(hfmc_array1) #单个预报分类ts评分</code></pre>
<pre><code>array([0.2 , 0.16666667, 0.5 ])</code></pre>
<pre><code class="language-python">mem.ts_hfmc(hfmc_grade1) #单个预报分级ts评分</code></pre>
<pre><code>array([1. , 0.25, 0.5 ])</code></pre>
<pre><code class="language-python">mem.ts_hfmc(hfmc_array2) #多个预报分类ts评分</code></pre>
<pre><code>array([[0.21052632, 0.18181818, 0.41935484],
[0.25 , 0.13333333, 0.40740741],
[0.08 , 0.06666667, 0.23333333],
[0.125 , 0.13333333, 0.24137931],
[0.18181818, 0.14285714, 0.31034483]])</code></pre>
<pre><code class="language-python">mem.ts_hfmc(hfmc_grade2) #多个分类分级ts评分</code></pre>
<pre><code>array([[1. , 0.4 , 0.41935484],
[1. , 0.25 , 0.40740741],
[1. , 0.125 , 0.23333333],
[1. , 0.25 , 0.24137931],
[1. , 0.28571429, 0.31034483]])</code></pre>
<pre><code class="language-python">mem.ets_hfmc(hfmc_array1)</code></pre>
<pre><code>array([0.02587519, 0.08045977, 0.17647059])</code></pre>
<pre><code class="language-python">mem.ets_hfmc(hfmc_array2)</code></pre>
<pre><code>array([[ 0.04458599, 0.10891089, 0.05263158],
[ 0.07120743, 0.01515152, 0.11111111],
[-0.11650485, -0.04283054, -0.06976744],
[-0.07142857, 0.01515152, -0.04761905],
[-0.00699301, 0.03420523, 0. ]])</code></pre>
<pre><code class="language-python">mem.bias_hfmc(hfmc_array1)</code></pre>
<pre><code>array([1.18181818, 0.55555556, 1.1 ])</code></pre>
<pre><code class="language-python">mem.bias_hfmc(hfmc_array2)</code></pre>
<pre><code>array([[1.09090909, 0.44444444, 1.2 ],
[1.27272727, 0.88888889, 0.9 ],
[1.45454545, 0.77777778, 0.85 ],
[1.45454545, 0.88888889, 0.8 ],
[1.36363636, 0.77777778, 0.9 ]])</code></pre>
<pre><code class="language-python">mem.far_hfmc(hfmc_array1)</code></pre>
<pre><code>array([0.69230769, 0.6 , 0.36363636])</code></pre>
<pre><code class="language-python">mem.far_hfmc(hfmc_array2)</code></pre>
<pre><code>array([[0.66666667, 0.5 , 0.45833333],
[0.64285714, 0.75 , 0.38888889],
[0.875 , 0.85714286, 0.58823529],
[0.8125 , 0.75 , 0.5625 ],
[0.73333333, 0.71428571, 0.5 ]])</code></pre>
<pre><code class="language-python">mem.mr_hfmc(hfmc_array1)</code></pre>
<pre><code>array([0.63636364, 0.77777778, 0.3 ])</code></pre>
<pre><code class="language-python">mem.mr_hfmc(hfmc_array2)</code></pre>
<pre><code>array([[0.63636364, 0.77777778, 0.35 ],
[0.54545455, 0.77777778, 0.45 ],
[0.81818182, 0.88888889, 0.65 ],
[0.72727273, 0.77777778, 0.65 ],
[0.63636364, 0.77777778, 0.55 ]])</code></pre>
<p>上述基于中间结果的检验函数也可以处理更高维度的输入，在此不做进一步举例。 </p>
<h1>Murphy 评分</h1>
<p>Murphy 评分是一种用于三分类预报问题的评分方法，最常用的场景是中长期预报中关于 偏高、正常、偏低 的分类检验。</p>
<p>Murphy评分的基本思路是对预报观测列联表矩阵（3×3的矩阵）的每个元素进行权重求和。列联表中对角线的元素代表正确预报，应该是奖励，因此权重应该是正的，非对角线元素是错误预报，应该惩罚，因此权重应该是负的。权重系数实际上也是一个3×3的矩阵，Murphy的评分解决的关键问题是设计一套具有公平性和合理性的权重系数矩阵。关于murphy细论述可参考《预报检验-大气科学从业者指南》 4.3.2节</p>
<p>计算Murphy 评分步骤包括：</p>
<ol>
<li>根据实况和预报数据，以及等级划分阈值，统计列联表</li>
<li>根据3分类预报问题的实况气候概率p1,p2和p3，确定Murphy评分的评分权重矩阵（函数 get_s_array）</li>
<li>根据列联表和权重矩阵，计算Murphy 评分。</li>
</ol>
<h2>Murphy评分权重系数</h2>
<p><strong><font face="黑体" color=blue size = 5>get_s_array(p1,p2,p3,k1=-0.25,k2=-0.25)</font></strong><br />
根据3分类预报问题的实况气候概率p1,p2和p3， 确定Murphy评分的评分权重矩阵。
评分权重的确定原理简单叙述如下：<br />
3分类预报的完整检验信息可由一个3×3的列联表表示。列联表中的元素p_ij表示 观测为第i类,预报为第j类的样本数。
总的评分等于列联表每个元素×一个权重系数后求和得到 sum(p_ij * s_ij)。 为此需要确定9个评分权重系数。
一般来说对角线上的权重系数s_ii 应该是正的，表示奖励， 非对角线上的s_ij是负的表示惩罚。那这9个评分权重该取为多少，
用它们做评价才是公平的？ 为此主要考虑2点：</p>
<ol>
<li>一是如果全部报对（列联表的元素全部集中在对角线上），评分应为1，</li>
<li>随机预报的评分应该为0，这2个条件可以构成关于权重系数的4个方程构成的方程组，这个方程组中包含3个分类的实况概率p1、p2、p3，
构成的已知参数。如果考虑这个评分矩阵是对称的，则一共有6个待定系数，因此它是一个欠定问题。为此有两个评分权重系数必须事先指定，
为此令 k1 = s12, k2 = s23， 在此基础上可以计算出其它4个系数 s11,s22,s33,s13.</li>
<li>计算出的评分权重还需满足： s12 < s11,s12<s22, s23<s22,s23<s33, s13<s12,s13<s23。 根据Murphy的分析，k1和k2的取值范围为(-0.5,0.5)时这些不等式条件都能得到满足。</li>
</ol>
<table>
<thead>
<tr>
<th style="text-align: left;">参数</th>
<th style="text-align: left;">说明</th>
</tr>
</thead>
<tbody>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>p1</font></strong></td>
<td style="text-align: left;">实况第1分类的概率</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>p2</font></strong></td>
<td style="text-align: left;">实况第2分类的概率</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>p3</font></strong></td>
<td style="text-align: left;">实况第3分类的概率</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>k1</font></strong></td>
<td style="text-align: left;">权重系数s12</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>k2</font></strong></td>
<td style="text-align: left;">权重系数s23</td>
</tr>
<tr>
<td style="text-align: left;"><font face="黑体" color=blue size=5>return</font></td>
<td style="text-align: left;">3×3的评分权重矩阵</td>
</tr>
</tbody>
</table>
<p><strong> 应用示例</strong> </p>
<pre><code class="language-python">s_array = mem.get_s_array(0.2,0.5,0.3,k1 = -0.5,k2 = -0.25)
print(s_array)</code></pre>
<pre><code>[[ 2.6 -0.5 -0.9 ]
[-0.5 0.35 -0.25 ]
[-0.9 -0.25 1.01666667]]</code></pre>
<h2>给定权重下的评分</h2>
<p><strong><font face="黑体" color=blue size = 5>murphy_ctable(p_array,s_array)</font></strong> </p>
<table>
<thead>
<tr>
<th style="text-align: left;">参数</th>
<th style="text-align: left;">说明</th>
</tr>
</thead>
<tbody>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>ctable_array</font></strong></td>
<td style="text-align: left;">3×3的列联表</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>s_array</font></strong></td>
<td style="text-align: left;">3×3的权重系数矩阵</td>
</tr>
<tr>
<td style="text-align: left;"><font face="黑体" color=blue size=5>return</font></td>
<td style="text-align: left;">评分，实数</td>
</tr>
</tbody>
</table>
<p><strong> 应用示例</strong> </p>
<pre><code class="language-python">ob = np.random.randn(2,10000)
fo1 = ob + np.random.randn(2,10000) *1
ctable_array = mem.contingency_table_multicategory(ob,fo1,grade_list = [-0.5,0.5])
print(ctable_array) </code></pre>
<pre><code>[[ 4413. 2415. 415. 7243.]
[ 1409. 2754. 1351. 5514.]
[ 441. 2457. 4345. 7243.]
[ 6263. 7626. 6111. 20000.]]</code></pre>
<p>mem.contingency_table_multicategory返回的列联表还包含累加量，因此是4×4的。</p>
<pre><code class="language-python">p_array = ctable_array[:-1,:-1]/ctable_array[1,1] #计算Murphy时，只需取3×3的部分,并除以总样本数，得到概率性的列联表。
score = mem.murphy_ctable(p_array,s_array)
print(score) </code></pre>
<pre><code>4.800550714112806</code></pre>
<h2>Murphy评分综合计算</h2>
<p><strong><font face="黑体" color=blue size = 5>murphy_score(ob,fo,grade_list = None,p1 = None,p2 = None,p3 = None,k1 = -0.25,k2 = 0.25)</font></strong> </p>
<p>根据原始输入数据，等级阈值，</p>
<table>
<thead>
<tr>
<th style="text-align: left;">参数</th>
<th style="text-align: left;">说明</th>
</tr>
</thead>
<tbody>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>ob</font></strong></td>
<td style="text-align: left;">实况数据, 任意维numpy数组</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>fo</font></strong></td>
<td style="text-align: left;">fo比Ob.shape多一维或者保持一致，fo.shape低维与ob.shape保持一致</td>
</tr>
<tr>
<td style="text-align: left;"><strong>grade_list</strong></td>
<td style="text-align: left;">如果该参数为None，观测或预报值出现过的值都作为分类标记，如果该参数不为None，它必须是一个从小到大排列的实数，以其中列出的数值划分出的多个区间作为分类标签。对于预报和观测值不为整数的情况，grade_list 不能设置为None</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>p1</font></strong></td>
<td style="text-align: left;">实况第1分类的概率</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>p2</font></strong></td>
<td style="text-align: left;">实况第2分类的概率</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>p3</font></strong></td>
<td style="text-align: left;">实况第3分类的概率</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>k1</font></strong></td>
<td style="text-align: left;">权重系数s12</td>
</tr>
<tr>
<td style="text-align: left;"><strong><font face="黑体" color=blue size = 5>k2</font></strong></td>
<td style="text-align: left;">权重系数s23</td>
</tr>
<tr>
<td style="text-align: left;"><font face="黑体" color=blue size=5>return</font></td>
<td style="text-align: left;">如果Fo和Ob的shape一致（即只有一种预报），返回结果为实数；如果Fo比Ob高出一维，则返回1维数组</td>
</tr>
</tbody>
</table>
<p><strong> 应用示例</strong> </p>
<pre><code class="language-python">ob = np.random.randn(2,10000)
fo2= ob + np.random.randn(3,2,10000) *1
score = mem.murphy_score(ob,fo2,grade_list = [-0.5,0.5]) # 根据数据自动计算p1,p2,p3，并且确定权重系数矩阵，再计算评分
print(score) </code></pre>
<pre><code>[1.80834036 1.78610498 1.7945716 ]</code></pre>
<pre><code class="language-python">#如果ob 和fo是离散的分类标记
n = 20
type_ = np.array([-1,0,1])
index_ob = np.random.randint(0, 3, size=n)
ob = type_[index_ob]
print(ob)</code></pre>
<pre><code>[-1 -1 0 -1 1 1 0 0 0 1 -1 -1 0 -1 1 -1 -1 0 -1 0]</code></pre>
<pre><code class="language-python">index_fo = np.random.randint(0, 3, size=n)
fo = type_[index_fo]
print(fo)</code></pre>
<pre><code>[-1 1 0 1 1 1 0 0 0 1 0 -1 1 0 0 1 -1 1 -1 -1]</code></pre>
<pre><code class="language-python">score = mem.murphy_score(ob,fo) # 根据数据自动计算p1,p2,p3，并且确定权重系数矩阵，再计算评分
print(score) </code></pre>
<pre><code>1.3241158642944355</code></pre>