How to calculate Probability of detection (hit rate) (POD),False alarm ratio (FAR), Accuracy (ACC) between two precipitation products.

27 次查看（过去 30 天）
shanka sharma
shanka sharma 2019-1-3
回答： shankar sharma 2022-8-18
i want to calculate POD, FAR, and ACC, between 387 station daily data sets with GSMAP precipitaion with similar dimension data sets, code and data sets are attached; data also contain some NaN values.
clc
clear all
load('pre_product.mat')
load('pre_gauge.mat')
for m=1:length(pre_gauge);%%(number of staitons)
aa=0;
bb=0;
cc=0;
dd=0;
pre_g_hit=[];
pre_p_hit=[];
pre_g_false=[];
pre_p_false=[];
pre_g_miss=[];
pre_p_miss=[];
for k=1:365;%%(day)
if (pre_gauge(k)>=0.1)&&(pre_product(k)>=0.1)
aa=aa+1;
pre_g_hit=[pre_g_hit;pre_gauge(k)];
pre_p_hit=[pre_p_hit;pre_product(k)];
elseif (pre_gauge(k)<0.1)&&(pre_product(k)>=0.1)
bb=bb+1;
pre_g_false=[pre_g_false;pre_gauge(k)];
pre_p_false=[pre_p_false;pre_product(k)];
elseif (pre_gauge(k)>=0.1)&&(pre_product(k)<0.1)
cc=cc+1;
pre_g_miss=[pre_g_miss;pre_gauge(k)];
pre_p_miss=[pre_p_miss;pre_product(k)];
elseif (pre_gauge(k)<0.1)&&(pre_product(k)<0.1)
dd=dd+1;
end
end
r0=corrcoef(pre_product,pre_gauge);
R(m)=r0(1,2);
bias(m)=sum(pre_product-pre_gauge)/sum(pre_gauge);
bias_hit(m)=sum(pre_p_hit-pre_g_hit)/sum(pre_gauge);
bias_false(m)=sum(pre_p_false-pre_g_false)/sum(pre_gauge);
bias_miss(m)=sum(pre_p_miss-pre_g_miss)/sum(pre_gauge);
RMSE(m)=(sum((pre_product-pre_gauge).^2)/length(pre_gauge))^0.5;
POD(m)=aa/(aa+cc);
FAR(m)=bb/(aa+bb);
ACC(m)=(aa+dd)/(aa+bb+cc+dd);
end
  3 个评论
显示 1更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

回答（1 个）

shankar sharma
shankar sharma 2022-8-18
Recently, i have solved the question,
here are the codes
pre_gauge=reshape(obs,387,1095); %% making station data timeseries (387 stations and 1095 days (3year))
pre_product=reshape(gpm,387,1095); %% making satellite data timeseries (387 stations and 1095 days (3year))
%%loop for each stations (it calculates for each station, a total of 387
%%staitons)
for m=1:387;%%(number of staitons)
aa=0;
bb=0;
cc=0;
dd=0;
pre_g_hit=[];
pre_p_hit=[];
pre_g_false=[];
pre_p_false=[];
pre_g_miss=[];
pre_p_miss=[];
for k=1:1095; %%(number of day)
if (pre_gauge(m,k)>=1)&&(pre_product(m,k)>=1)
aa=aa+1;
pre_g_hit=[pre_g_hit;pre_gauge(m,k)];
pre_p_hit=[pre_p_hit;pre_product(m,k)];
elseif (pre_gauge(m,k)<1)&&(pre_product(m,k)>=1)
bb=bb+1;
pre_g_false=[pre_g_false;pre_gauge(m,k)];
pre_p_false=[pre_p_false;pre_product(m,k)];
elseif (pre_gauge(m,k)>=1)&&(pre_product(m,k)<1)
cc=cc+1;
pre_g_miss=[pre_g_miss;pre_gauge(m,k)];
pre_p_miss=[pre_p_miss;pre_product(m,k)];
elseif (pre_gauge(m,k)<1)&&(pre_product(m,k)<1)
dd=dd+1;
end
end
%% formula for pod, far, acc
POD(m)=aa/(aa+cc);
FAR(m)=bb/(aa+bb);
ACC(m)=(aa+dd)/(aa+bb+cc+dd);
CSI(m)=aa/(aa+bb+cc);
FBI(m)=(bb+aa)/(aa+cc);
end
%% average for whole stations
ALL_POD=AA/(AA+CC)
ALL_FAR=BB/(AA+BB)
ALL_ACC=(AA+DD)/(AA+BB+CC+DD)
ALL_FBI=(AA+BB)/(AA+CC)
ALL_CSI=AA/(AA+BB+CC)

类别

在 Help Center 和 File Exchange 中查找有关 Applications 的更多信息

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by