%% Clear the WorkSpace and CommandWindow
clear all
close all
clc

%% read data from in situ pH
% data = [year BotDep Lon Lat month Depth T S DO N P Si insitupHT ANNpHT error];
data15 = xlsread('F:\paper result\pHT\Final result2\comparison volumn\six stations_water column.xlsx','15A');
data16 = xlsread('F:\paper result\pHT\Final result2\comparison volumn\six stations_water column.xlsx','16A');
data67 = xlsread('F:\paper result\pHT\Final result2\comparison volumn\six stations_water column.xlsx','67A');
data69 = xlsread('F:\paper result\pHT\Final result2\comparison volumn\six stations_water column.xlsx','69A');
data75 = xlsread('F:\paper result\pHT\Final result2\comparison volumn\six stations_water column.xlsx','75A');
data85 = xlsread('F:\paper result\pHT\Final result2\comparison volumn\six stations_water column.xlsx','85A');

data = [data15;data16;data67;data69;data75;data85];
targets = data(:,13); % insitu pH 
outputs = data(:,14); % ANN pH
X =[ones(size(outputs,1),1),targets];
MAE     = sum(abs(targets-outputs))/size(outputs,1);
RMSE    = sqrt( sum((targets-outputs).^2)/size(outputs,1) );
R2      = 1-sum((targets-outputs).^2)/sum((targets-mean(targets)).^2);
[b,bint,r,rint,stats] = regress(outputs,X);
forcast =  b(1) + b(2) * targets;
x = [7.7:0.1:8.4];
y = x;
figure;
scatter(data15(:,13),data15(:,14),'ro');
hold on;
scatter(data16(:,13),data16(:,14),'bo');
hold on;
scatter(data67(:,13),data67(:,14),'go');
hold on;
scatter(data69(:,13),data69(:,14),'rs');
hold on;
scatter(data75(:,13),data75(:,14),'bs');
hold on;
scatter(data85(:,13),data85(:,14),'gs');
axis([7.7,8.4,7.7,8.4]);%设置坐标轴上下限
set(gca,'FontName','Times New Roman','FontSize',12);%坐标轴的字体、大小
% title('Station A1-5','FontName','Times New Roman','FontSize',10);
legend('A1-5','A1-6','A6-7','A6-9','A7-5','A8-5','FontName','Times New Roman','FontSize',12);
xlabel('Observation pH_T','FontName','Times New Roman','FontSize',12)
ylabel('Retrieved pH_T from FVCOM output','FontName','Times New Roman','FontSize',12)
hold on;
plot(x,y,'k--');
text(7.75,8.35,'N = 84','FontName','Times New Roman','FontSize',12);
text(7.75,8.32,'MAE = 0.04','FontName','Times New Roman','FontSize',12);
text(7.75,8.29,'RMSE = 0.05','FontName','Times New Roman','FontSize',12);
text(7.75,8.26,'R^{2} = 0.71','FontName','Times New Roman','FontSize',12);
    
