-
Notifications
You must be signed in to change notification settings - Fork 2
/
GetNonParaFit.m
153 lines (121 loc) · 4.11 KB
/
GetNonParaFit.m
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
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
function [Cancel_Flag, Abunds, EMs, Fit_Quality] = GetNonParaFit(Xin, GS, Fit_Params)
%
% Function that returns the non-parametric endmember fit
%
% Input:
% Xin - nData x nVar matrix containing the observed data
% Fit_Params - 5 x 1 cell of parameters specifying fit options
%
% Output:
% Cancel_Flag - a flag to indicate whether (= 1) or not (= 0) the
% user cancelled the process
% Abunds - nData x nEnd matrix of end member abundances
% EMs - nEnd x nVar matrix of end members (NB not multiplied by A)
% Dist_Params - nEnd x nParams matrix of the best-fit distribution
% parameters
% Fit_Quality - cell array containing measures of the goodness of
% fit for the best-fit endmembers{[R2(Selected_EM,
% Mean_Angle], [Spec_R2, Spec_Angle]};
%
%% Process the inputs
X = Xin; % Assign the original data to X
[nData, nVar] = size(X);
% X(X<eps) = eps; % round off very small numbers and zeros (prevents errors)
% X = X./(sum(X,2)*ones(1,nVar)); % set sum of each row to 1
% The min and max End Members to fit
EM_Min = Fit_Params{1};
EM_Max = Fit_Params{2};
%% Loop through the endmembers
% Pre-allocate variables for speed
R2 = NaN(EM_Max,1);
Spec_R2 = NaN(nData, EM_Max);
Min_Spec_R2 = NaN(EM_Max, 1);
EM_R2 = NaN(EM_Max,1);
Spec_Angle = NaN(nData, EM_Max);
DataSet_Angle = NaN(EM_Max, 1);
% Convexity_Error = NaN(EM_Max, 1);
Stored_Xprime = cell(EM_Max, 1); % for storing the abundance
Stored_EMs = cell(EM_Max, 1); % for storing the end members
Stored_Abunds = cell(EM_Max, 1); % for storing the abundance
% tic
% Set up a waitbar to count the loop
h = waitbar(0,'Initializing....', 'Name', 'Calculating end member fits...',...
'CreateCancelBtn', 'setappdata(gcbf,''Cancelled'',1)');
setappdata(h,'Cancelled',0)
Cancel_Flag = 0;
for k = EM_Min:EM_Max
% Check for Cancel button press
if getappdata(h,'Cancelled')
Cancel_Flag = 1;
break;
end
% Update the waitbar and continue
waitbar((k-1)/(EM_Max), h, strcat('Fitting ', sprintf(' %d', k), ' end members....'))
tmp_EM = HALS_NMF(X, k, 5e3, 10, [5, 0], 0);
% Check convexity and adjust b2 if needed
[tmp_EM, tmp_Abunds, Xprime, Convexity] = Find_b2(X, tmp_EM);
% Sort the EMs
[tmp_EM, Sinds] = sortEMs(tmp_EM, GS, 'Median');
tmp_Abunds = tmp_Abunds(:,Sinds);
Stored_Xprime(k) = {Xprime};
Stored_EMs(k) = {tmp_EM};
Stored_Abunds(k) = {tmp_Abunds};
%% Get the correlations and update waitbar
R2(k) = GetR2(Xprime(:), X(:));
Spec_R2(:, k) = GetR2(X', Xprime')';
Min_Spec_R2(k) = min( Spec_R2(:, k) );
DataSet_Angle(k) = GetAngles(X(:), Xprime(:));
Spec_Angle(:,k) = GetAngles(X, Xprime);
% Convexity_Error(k) = GetConvexityError(X, tmp_EM);
if k >1
r = GetR2(tmp_EM');
r = r - diag(diag(r));
EM_R2(k) = max(max(r.^2));
end
% Check for cancelling
if getappdata(h,'Cancelled')
Cancel_Flag = 1;
break;
end
end
delete(h)
% toc
if Cancel_Flag == 1
% User has cancelled so return emptiness
Abunds=[];
EMs=[];
Fit_Quality=[];
return
end
% Get the selected number of end members
Transfer.DataSet_R2 = R2;
Transfer.Spec_R2 = Spec_R2;
Transfer.DataSet_Angle = DataSet_Angle;
Transfer.Spec_Angle = Spec_Angle;
Transfer.EM_Min = EM_Min;
Transfer.EM_Max = EM_Max;
Transfer.EM_R2 = EM_R2;
% Transfer.Convexity = Convexity_Error;
Return = Select_EndMembers('DataTransfer', Transfer);
Cancel_Flag = Return.Cancel_Flag;
if Cancel_Flag == 1
% User has cancelled so return emptiness
Abunds=[];
EMs=[];
Fit_Quality=[];
return
end
Selected_EM = Return.EM;
Xprime = Stored_Xprime{Selected_EM};
EMs = Stored_EMs{Selected_EM};
Abunds = Stored_Abunds{Selected_EM};
% Recalculate the fit quality stats
R2 = GetR2(Xprime(:), X(:));
Spec_R2 = GetR2(X', Xprime')';
DataSet_Angle = GetAngles(X(:), Xprime(:));
Spec_Angle = GetAngles(X, Xprime);
r = GetR2(EMs');
r = r - diag(diag(r));
EM_R2 = max(max(r.^2));
Fit_Quality = {[R2, DataSet_Angle, EM_R2], [Spec_R2, Spec_Angle]};
% toc