% rename ps_monod.m 10/18/22
% 2018_PS3_monod_fit.m, mak saito 3/4/18
% uses monod.m function and nnleasqr.m and dfdo.m for fitting program
%  Modified from CoPhaeoMM4.m (Feb16, 2006) 3/4/18
%  Phaeocystis Culture Data - Simple Experiment
%  blank corrected no carryover

clear
clf

% Culture data tables for analysis

% data from Transfer 2 averaged of replicates A and B 
% Co2+  Zn2+  Avgu  

% add columns of Sunda and Huntsman Data here and remove % symbol to make code active
%Codata=[];
%Zndata=[];

% Phaeocystis data - delete or add % to beginning of each line to deactivate
% Co2+  Zn2+  Avgu  
Codata=[
-13.19	-13.36	0.04
-12.86	-13.36	0.05
-12.51	-13.36	0.07
-12.10	-13.36	0.09
-11.60	-13.36	0.10
-11.13	-13.36	0.10
-10.62	-13.36	0.12
-10.15	-13.36	0.14];
% data from T2
% Co2+  Zn2+  Avgu  
Zndata=[
-13.19	-13.12	0.04
-13.19	-12.82	0.04
-13.19	-12.44	0.06
-13.19	-11.96	0.13
-13.19	-11.49	0.20
-13.19	-10.98	0.20
-13.19	-10.51	0.26];


co=Codata(:,1);
con=10.^(co);
uco=Codata(:,3);


zn=Zndata(:,2);
znn=10.^(zn);
uzn=Zndata(:,3);

figure(1)
clf

% PLOT ORIGINAL DATASET 
subplot(121) 
semilogx(znn,uzn,'*r-');
hold on
semilogx(con, uco,'+b-');
title('original data')
ylabel('Chl specific growth rate (per day)');
xlabel('Zinc or Cobalt concentration (nM)');

% Non-least Squares analysis Zn data
S = znn;
u = uzn;		
pin=[0.26 1e-13];  % initial guesses for umax and Km
% 	run non-linear least squares function for Zn
[f,p,kvg,iter,corp,covp,covr,stdresid,Z,r2]=nlleasqr(S,u,pin,'monod');

%	plot original data and fit
subplot(122) 
C=((1:10000).*1e-14)';
y=monod(C,p);
semilogx(S,u,'+'); hold on
semilogx(C,y);
xlabel('Zinc concentration (nM)');
ylabel('Chl specific growth rate (per day)');
%axis([1e-15 1e-9 0 .3]);
zn_umax=p(1)     	% 	record the solved umax	
zn_K=p(2)		%	the solved Ks
error=diag(covp).^(0.5)	% the parameter errors are the square roots of the diagonal
								% of the covp matrix
r2								% r2 error value


% Non-least Squares analysis Co data
S = con;
u=uco;
pin=[0.20 1e-13];  % initial guesses for umax and Km
% 	run non-linear least squares function
[f,p,kvg,iter,corp,covp,covr,stdresid,Z,r2]=nlleasqr(S,u,pin,'monod');
C=((1:10000).*1e-14)';
y=monod(C,p);
semilogx(S,u,'*'); hold on
semilogx(C,y);
xlabel('Zinc or Cobalt concentration (nM)');
ylabel('Chl specific growth rate (per day)');

%axis([1e-15 1e-9 0 .3]);
co_umax=p(1)   	% 	 the solved umax
co_Km=p(2)       %    the solved Ks
error=diag(covp).^(0.5)	% the parameter errors are the square roots of the diagonal
								% of the covp matrix
r2								% r2 error value