-
Notifications
You must be signed in to change notification settings - Fork 0
/
MAhpcDrift.asv
76 lines (76 loc) · 2.11 KB
/
MAhpcDrift.asv
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
function [f,g]=MAhpcDrift%(lscan,fluscapeLocations)
thismany=1;%Random ICs in fSMR - loop length
tauend=1000;
burn=500;
years=tauend-burn;
isdual=0;
solvetype=2;
numseed=10^(-5);
eps=(.3:.1:.7);
leps=length(eps);
%
load('forMAhpc.mat')
[gamma,NN,n,nbar,na,NNbar,NNrep,minNind,maxNind,maxN,Kbar,K1,Cbar,betaS,betaI,betaD,beta3,ages0]=prepFluAgeLocs(C,Qeven,0,1);
%}
%{
load('fluscape.mat')
[Df,xf,yf,Ds,Is]=fluscapeOnlyLocalAverage(lscan,fluscapeLocations);
fscapeind=find(Is);
[lscanNew,r]=fluscapeNNr(Df,Ds);
[gamma,NN,n,nbar,na,NNbar,NNrep,minNind,maxNind,maxN,Kbar,K1,Cbar,betaS,betaI,betaD,beta3,ages0]=prepFluAgeLocsFscape(lscanNew(fscapeind),r(fscapeind,fscapeind),0,1);
%}
[fp,gp,~]=finalSizeMulti(gamma,n,nbar,na,NN,NNbar,NNrep,minNind,maxNind,maxN,Kbar,K1,Cbar,betaS,betaI,betaD,beta3,isdual,solvetype,numseed,0,0,1);
%
%X=zeros(n,leps,thismany);
%Y=X;
X=nan(leps,years);
Y=X;
thresh=.005;
yy=1:years;
for i=1:leps
epsi=eps(i);
for j=1:thismany
[f,g,~]=finalSizeMulti(gamma,n,nbar,na,NN,NNbar,NNrep,minNind,maxNind,maxN,Kbar,K1,Cbar,betaS,betaI,betaD,beta3,isdual,solvetype,numseed,epsi,1,500);
%
fx=f(:,burn+1:end); %g=g(burn+1:end);
fx1=fx;
%fsum=sum(fx,1);
fsum=max(fx,[],1);
fx(:,fsum<thresh)=[];
yyi1=yyl; yyi1(fsum<thresh)=[];
lf=size(fx,2);
%
gx=g(:,burn+1:end); %g=g(burn+1:end);
gx1=gx;
%gsum=sum(gx,1);
gsum=max(gx,[],1);
gx(:,gsum<thresh)=[];
yyi2=yy; yyi2(fsum<thresh)=[];
lg=size(gx,2);
%
XX=zeros(n,lf);
YY=zeros(n,lg);
for k=1:lf
%{
xi=[(1:n)',fx(:,k)];
xi=sortrows(xi,2);
XX(:,k)=xi(:,1);
%}
cck=corrcoef(fx(:,k),fp);
X(i,k)=
end
for k=1:lg
%{
yi=[(1:n)',gx(:,k)];
yi=sortrows(yi,2);
YY(:,k)=yi(:,1);
%}
end
X(:,i,j)=var(XX,0,2);
Y(:,i,j)=var(YY,0,2);
end
end
f=X;
g=Y;
save('MA1drift','X','Y')
end