(Die Zahlenwerte für g und h sind nur ein Beispiel)
Fehlermeldung:
Error using fit>iFit (line 415)
NaN computed by model function, fitting cannot continue.
Try using or tightening upper and lower bounds on coefficients.
Error in fit (line 109)
[fitobj, goodness, output, convmsg] = iFit( xdatain, ydatain, fittypeobj, ...
Error in testfitR2012a (line 9)
gfit = fit(g',h',f1);
Mit der 2009b Version von Matlab funktioniert dieses Programm jedoch fehlerfrei mit folgendem Ergebnis:
gfit =
General model:
gfit(x) = a*(1+c*x^2)*exp(-((x^2)/(2*s^2)))
Coefficients (with 95% confidence bounds):
a = 3.249 (2.715, 3.784)
c = 8.19e-008 (fixed at bound)
s = 8.592 (-6.586, 23.77)
Wie bekomme ich nun mein altes Programm auf der neuen Matlabversion zum laufen?
Hallo Harald,
ich habe ein ähnliches Problem wie der Herr damals.
Es ist nicht mein Code, sondern der meines Betreuers, mit allen vorigen Daten hat es aber wunderbar geklappt.
Code:
[YY=Meanmeanclusterhisty;
XXX=Meanmeanclusterhistx;
YY = YY./max(YY);
ff=fittype('gauss1');
fopts=fitoptions('gauss1');
c = find(YY~=0);
xx = XXX(c);
yy = YY(c);
[cfun,gof,output]...
=fit(xx,yy,ff);
=fit(Meanmeanclusterhistx,Meanmeanclusterhisty,ff);
figure;
hold on
plot(XXX,YY)]
Die Verteilung, siehe Bild im Anhang, sieht soweit eigentlich gut aus. Scheinbar habe auch ich ein Problem mit einem negativen s bzw c. Folgende Fehlermeldung:
Error using fit>iFit (line 340)
Inf computed by model function,
fitting cannot continue.
Try using or tightening upper and
lower bounds on coefficients.
Error in fit (line 108)
[fitobj, goodness, output, convmsg]
= iFit( xdatain, ydatain,
fittypeobj, ...
Error in min_max_frequ_mean_v2 (line
811)
[cfun,gof,output]...
Kannst du mir weiterhelfen und verraten, wie ich das Problem gelöst bekomme?
vielen Dank!
Hallo Harald,
entschuldige, ich werde die Daten anhängen...
Hier der ganze Code (sorry, ziemlich lang) Ich hänge ihn auch nochmal als mat.file an (min_max_frequ_mean_v2)
Code:
clearall
kkk1=0;
aa=0;
while aa==0 clear FFinv FFinv1 clustermitos1 Ensemble signal3 FFmean FFstd Ensemble_peak peakcutoff cmmi cmma
clear data data2 Ensemble signal2 totpeaks
clear muf_rw ZZx ZZy L data data2 Maxhist dminus dplus
kkk1=kkk1+1;
prompt={'Proceed with cell' };
defans={num2str(kkk1)};
fields = {'one'};
options.Resize='on';
options.WindowStyle='normal';
info = inputdlg(prompt, 'Go to', 1, defans,options);
if ~isempty(info) %see if user hit cancel
info = cell2struct(info,fields);
[filename,pathname]= uigetfile(['C:\CVRC\MitochondrialMetabolism\*.mat'],'Choose the _frequ_dist file:');
ifexist([pathname filename])==2% check whether there is a file
fid=fopen([pathname filename]); % open file
position=0;
while ~feof(fid);
data=load([pathname filename],'params','mitosincluster');
file{kkk1}=filename(1:6);
fclose(fid);break;
end;
else display('You suck, take an existing file!') end
rw=data.params.rw;
fs=data.params.fs;
bb{kkk1}=data.params.bb;
bbend{kkk1}=data.params.bbend;
mitosincluster{kkk1}=data.mitosincluster;
clear data
prompt={'Name following file: '};
defans={filename(1:10)};
fields = {'name'};
options.Resize='on';
options.WindowStyle='normal';
info = inputdlg(prompt, 'Mito#: ', 1, defans,options);
if ~isempty(info) %see if user hit cancel
info = cell2struct(info,fields);
filename3 = info.name;
end
k1=1;
[filename1,pathname1]= uigetfile(pathname,'Choose the _wavefrequs_final:');
ifexist([pathname1 filename1])==2% check whether there is a file
fid=fopen([pathname1 filename1]); % open file
position=0;
while ~feof(fid);
data=load([pathname1 filename1]);
fclose(fid);break;
end;
else display('You suck, take an existing file!') end [filename1,pathname1]= uigetfile(pathname,'Choose the _mat:');
ifexist([pathname1 filename1])==2% check whether there is a file
fid=fopen([pathname1 filename1]); % open file
position=0;
while ~feof(fid);
data2=load([pathname1 filename1]);
fclose(fid);break;
end;
else display('You suck, take an existing file!') end
clear ZZ
prompt={'Sampling frequency [sec]','Auto CC cutoff for maxpeak','rwsize','Spacing','Peak cutoff','CC cutoff','Alternative CC cutoff', 'Alternative2 CC cutoff','Intervall around peak [mHz]:','Peak-Intervall cutoff (innercutoff):','Normalized Mito cutoff: '};
defans={num2str((1./fs)*1000),num2str(0.5),num2str(rw),num2str(0.1),num2str(0.05),num2str(0.95),num2str(0.95),num2str(0.95),num2str(0.5),num2str(0.5),num2str(0.05)};
fields = {'fs','aCC','rw','name','name2','name3','name4','name5','roundseg','innercutoff','mitocutoff'};
options.Resize='on';
options.WindowStyle='normal';
info = inputdlg(prompt, 'Define parameters: ', 1,defans, options);
if ~isempty(info) %see if user hit cancel
info = cell2struct(info,fields);
seg=str2num(info.name);
peakcutoff1=str2num(info.name2);
ccut=str2num(info.name3);
cccut=str2num(info.name4);
cccut_alt2=str2num(info.name5);
rw=str2num(info.rw);
auto_CC_cutoff=str2num(info.aCC);
fs=1/str2num(info.fs).*1000; %%CAVE: fs in [mHz]!!
roundseg=str2num(info.roundseg);
innercutoff=str2num(info.innercutoff);
mito_cutoff=str2num(info.mitocutoff);
end [IND,num]=size(data.FF);
for k1=1:num
prompt={'Frames to initiation: ','Stop-frames: ','Upper-frequ-cutoff: ','Lower-frequ-cutoff: '};
defans={num2str(bb{kkk1}{k1}),num2str(bbend{kkk1}{k1}),... num2str(round2(fs./data.params.upper_cutoff{k1},log10(1/data.params.interpol_seg))),... num2str(round2(fs./(data.params.lower_cutoff{k1}),log10(1/data.params.interpol_seg)))};
fields = {'name','name2','name3','name4'};
options.Resize='on';
options.WindowStyle='normal';
info = inputdlg(prompt, 'Frames: ', 1, defans,options);
if ~isempty(info) %see if user hit cancel
info = cell2struct(info,fields);
bb{kkk1}{k1} = str2num(info.name);
bbend{kkk1}{k1}=str2num(info.name2);
frequUcut{k1}=str2num(info.name3);
frequLcut{k1}=str2num(info.name4);
end end
% seg=0.05; %take spacing as scale resolution dj (see wavelet) % peakcutoff=0.05 %%cut off of bigger peaks
for k1=1:num
ll=frequLcut{k1}:data.params.interpol_seg:frequUcut{k1};
vec=round(ll./seg);
for a=1:length(data.FF{2,1})
a
if ccc<vec(end)-seg %%make algorithm faster by not taking all frequ-scale points
ZZx{a,k1}=vec(1:find(vec==ccc+1)).*seg;
else
ZZx{a,k1}=vec.*seg;
end
ZZZx{a,k1,kkk1}=ZZx{a,k1};
% ZZx{a,k1}=interp1(x,y,xx); %%interpolated x-frequency-scale
ZZy{a,k1}=hist(muf_rw{a,k1,1},ZZx{a,k1}); %%number of corresponding mitos
ZZZy{a,k1,kkk1}=ZZy{a,k1};
end end
% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % %%%standard deviation and mean %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for k1=1:num
for i=1:length(data.FF{2,1})
c=find(FFinv{i,k1}>0);
FFinv1{i,k1}=FFinv{i,k1}(c);
cc=find(FFinv1{i,k1}>ZZx{i,k1}(1) & FFinv1{i,k1}<ZZx{i,k1}(end));
FFmean{k1}(i)=mean(FFinv1{i,k1}(cc));
FFstd{k1}(i)=std(FFinv1{i,k1});
end end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% %%now classify separated peaks and order them: %%%%%%%%%%%%%%%%%%%%%%%%%%%
for k1=1:num
for a=1:length(data.FF{2,1})
a
% peakcutoff{a,k1}=0.1.*; for kk=1:length(ZZx{a,k1}) %%CAVE: cut-off still in scale
Ensemble.bars{a,kk,k1}=... find(muf_rw{a,k1}>=ZZx{a,k1}(kk)-seg./2 & ...
muf_rw{a,k1}<=ZZx{a,k1}(kk)+seg./2); %%giving mitos in each barfor each RW-frame
end
Maxhist{a,k1}=min(find(ZZy{a,k1}(2:end-1)==max(ZZy{a,k1}(2:end-1))))+1; %find maximum of ZZy frequ distribution
% peakcutoff{a,k1}=FFmean{k1}(a)+FFstd{k1}(a);
peakcutoff{a,k1}=0.05.*ZZy{a,k1}(Maxhist{a,k1});
% clear p ZZy_n % p{k1}(1)=1; %set correlation coeff 0 % p{k1}(2)=0; % ZZy_n{a,k1}=ZZy{a,k1}; % Maxhistt2=min(find(ZZy_n{a,k1}(2:end-1)==min(ZZy_n{a,k1}(2:end-1))))+1; %find minimal point in Distr. as start for while loop % % % if a<rw/2 % Maxhist{a,k1}=min(find(ZZy_n{a,k1}(2:end)==max(ZZy_n{a,k1}(2:end-1))))+1; % % elseif a>=rw/2 & a<=length(data.gpmean{2,1})-rw/2 % z=0; % for k=1:2 % if p{k1}(2)<=auto_CC_cutoff %%start looking for first peak, that has Auto-CC > auto_CC_cutoff (0.75) % p{k1}(1)=p{k1}(2); % z=z+1; % ZZy_n{a,k1}(Maxhistt2)=0; % Maxhistt=Maxhistt2; % Maxhistt2=min(find(ZZy_n{a,k1}(2:end-1)==max(ZZy_n{a,k1}(2:end-1))))+1; %%determine max-frequ % if length(Ensemble.bars{a,Maxhistt2,k1})<=innercutoff.*length(Ensemble.bars{a,Maxhistt,k1}) % p{k1}(2)=0; % else % c=Ensemble.bars{a,Maxhistt2,k1}; % ll=(max(a-rw./2+1,1):min(length(data.gpmean{2,1}),a+rw./2)); %set window rw for inner CC % Z{a}=zeros(1,length(ll)); % for kk=1:length(c) % Z{a}=Z{a}+data.gpmean{c(kk)+1,k1}(ll); %%take mean of tmre signal of % end %% tmRe signal, the "+1" stems from the fact that the first entry in gpmean is 0 (compare wavelet_analysis_v2) % Z{a}=Z{a}./length(c); % LUz=median(Z{a}); % R{a,k1}=zeros(2,2); % for kk=1:length(c) % LUm=median(data.gpmean{c(kk)+1,k1}(ll)); % CORRCOEF{kk,k1}=... %%do CC % corrcoef(Z{a}-LUz,data.gpmean{c(kk)+1,k1}(ll)-LUm);%%get corr coeff and check in while loop if it satisfies the auto_CC_cutoff % R{a,k1}=R{a,k1}+CORRCOEF{kk,k1}; % end % R{a,k1}=R{a,k1}./length(c); % p{k1}(2)=R{a,k1}(2,1); % end % end % end % % if z==2 & p{k1}(2)<=p{k1}(1) % Maxhist{a,k1}=Maxhistt; % elseif z==2 & p{k1}(2)>p{k1}(1) & ... % length(Ensemble.bars{a,Maxhistt,k1})>=innercutoff.*length(Ensemble.bars{a,Maxhistt,k1}) % Maxhist{a,k1}=Maxhistt2; % else % Maxhist{a,k1}=Maxhistt2; % end % % elseif a>length(data.gpmean{2,1})-rw/2 % Maxhist{a,k1}=min(find(ZZy_n{a,k1}(2:end-1)==max(ZZy_n{a,k1}(2:end-1))))+1; % end
% for kk=max(Maxhist{a,k1}-(roundseg)./seg,2):1:min(Maxhist{a,k1}+(roundseg)./seg,length(ZZx{a,k1})-1) %%check in 1mHz interval around maxpeak if there are similar high peaks and take them in % if length(Ensemble.bars{a,kk,k1})>=innercutoff.*length(Ensemble.bars{a,Maxhist{a,k1},k1}) % test(kk-max(Maxhist{a,k1}-(roundseg)./seg,2)+1)=1; % else % test(kk-max(Maxhist{a,k1}-(roundseg)./seg,2)+1)=0; % end % end % ccmin=min(min(find(test>0))+max(Maxhist{a,k1}-(roundseg)./seg,2)-1,Maxhist{a,k1}); % ccmax=max(max(find(test>0))+max(Maxhist{a,k1}-(roundseg)./seg,2)-1,Maxhist{a,k1});
for k1=1:num %%determine the mean signal Ensemble.signalfor each max peak
for a=1:length(data.FF{2,1})
a
ll=(max(a-rw./2+1,1):min(length(data.gpmean{2,1}),a+rw./2));
for k=1:Ensemble.peak_no{a,k1}
c=find(Ensemble.peak_L{a,k1}==k); %%for each coherent cluster peak find corresponding histogram bars
Ensemble.signal{a,k1,k}=zeros(1,length(ll));
z=0;
for kk=1:length(c)
cc=Ensemble.bars{a,c(kk),k1}; %%find mitos in the bars for whole cluster peak
if Maxhist{a,k1}==c(kk)
Ensemble.maxpeak{a,k1}=k; %%find maximum cluster peak
end for kkk=1:length(cc)
z=z+1;
% Ensemble.signal{a,k1,k}=Ensemble.signal{a,k1,k}+smooth(data.gpmean{cc(kkk)+1}(a:a+127),80,'rloess')';
Ensemble.signal{a,k1,k}=Ensemble.signal{a,k1,k}+data.gpmean{cc(kkk)+1,k1}(ll); %%getmean signal forall cluster peaks end end
Ensemble.signal{a,k1,k}=Ensemble.signal{a,k1,k}./z; %divide by Number of Mitos in cluster peak
L{a,k1,k}=z; %%number of mitos in clusterpeak k
end
for k=1:length(ZZx{a,k1})
Ensemble.signal2{a,k1,k}=zeros(1,length(ll));
z=0;
for kk=1:length(Ensemble.bars{a,k,k1})
Ensemble.signal2{a,k1,k}=Ensemble.signal2{a,k1,k}+data.gpmean{Ensemble.bars{a,k,k1}(kk)+1,k1}(ll);
z=z+1;
end if z~=0
Ensemble.signal2{a,k1,k}=Ensemble.signal2{a,k1,k}./z;
end end
[maxtab,mintab]=peakdet(ZZy{a,k1},0.5); %%findmax and min of histogram
ss=interp1(maxtab(:,1),maxtab(:,2),1:length(ZZx{a,k1})); %%interpolate between max over range of histogram
Diff=diff(ss); %%find difference function
pCl=1;
MM=min(find(Ensemble.peak_L{a,k1}==Ensemble.maxpeak{a,k1}));
k=max(MM-2,2); %beginning of maximum cluster peak
MaP=Ensemble.maxpeak{a,k1};
while pCl>=ccut & k>1
k=max(MM-2,2);
while k>1 & Diff(k)>=0
k=k-1;
end
Minbreakleft=k; %%find first Min. of interpolated maximum function to the left
while Ensemble_peak{a,k1}(k)==0 & k>=2
k=k-1;
end %%find first peak after Min to the left that is above cutoff (i.e. belonging to a cluster peak)
eins=k;
k=Minbreakleft;
while k>1 & Diff(k)<=0 %%find first Max after Min to the left
k=k-1;
end if ZZy{a,k1}(k)>ZZy{a,k1}(k+1)
zwei=k;
else
zwei=k+1;
end
if eins==1 & zwei==1
Nextsignifpeakleft=2;
pCl=0;
Peaknol=1;
Minbreakleft=2;
else
Nextsignifpeakleft=max(eins,zwei); %%take Max. of left max and first peak above cutoff and check peak
% clear Peaknol
Peaknol=-1;
for k=1:Ensemble.peak_no{a,k1}
c=find(Ensemble.peak_L{a,k1}==k);
for kk=1:length(c) if c(kk)==Nextsignifpeakleft
Peaknol=k; %%find corresponding cluster peak
end end end
EE=zeros(1,length(ll));
zz=0;
if Peaknol~=-1 & Peaknol<MaP
for k=Peaknol+1:MaP
EE=EE+Ensemble.signal{a,k1,k}.*L{a,k1,k};
zz=zz+L{a,k1,k}; %%getmean signal of corresponding cluster peak
end else
EE=Ensemble.signal2{a,k1,Nextsignifpeakleft}./length(Ensemble.bars{a,Nextsignifpeakleft,k1});
end
if MaP>Peaknol & Peaknol~=-1
EE=EE./zz;
Sigmax=EE;
Sigleft=Ensemble.signal{a,k1,Peaknol};
LUmax=median(Sigmax);
LUl=median(Sigleft);
CORRCOEFl=corrcoef(Sigmax-LUmax,Sigleft-LUl);
pCl=CORRCOEFl(2,1); %%check if corresponding cluster peak is highly correlated to maximum cluster peak
MM=min(find(Ensemble.peak_L{a,k1}==Peaknol));
MaP=Peaknol;
elseif Peaknol==-1
EEmax=zeros(1,length(ll));
z=0;
for k=Minbreakleft:Maxhist{a,k1}
EEmax=EEmax+Ensemble.signal2{a,k1,k};
z=z+1;
end
EEmax=EEmax./z;
Sigmax=EEmax;
Sigleft=EE;
LUmax=median(Sigmax);
LUl=median(Sigleft);
CORRCOEFl=corrcoef(Sigmax-LUmax,Sigleft-LUl);
pCl=CORRCOEFl(2,1); %%check if corresponding cluster peak is highly correlated to maximum cluster peak
c=find(ZZy{a,k1}>0);
k=Nextsignifpeakleft-1;
while k>1 & ZZy{a,k1}(k)==0
k=k-1;
end
MM=k;
% MaP=Nextsignifpeakleft; else
pCl=0;
Peaknol=Peaknol-1;
end end end if Peaknol==1 & Minbreakleft<=min(find(Ensemble.peak_L{a,k1}==1)) %%determine start-bar in Histogram by also considering border conditions
cl=Minbreakleft;
elseif Peaknol==Ensemble.peak_no{a,k1}
cl=min(find(Ensemble.peak_L{a,k1}==Peaknol));
k=1;
while cl>=mintab(k,1)
k=k+1;
end
cl=max(2,mintab(k-1,1));
elseif Peaknol==-1
cl=Minbreakleft;
else
Peaknol=Peaknol+1;
cl=min(Minbreakleft,min(find(Ensemble.peak_L{a,k1}==Peaknol)));
k=1;
if cl>mintab(1,1) while cl>=mintab(k,1)
k=k+1;
end if k==1
cl=2;
else
cl=max(2,mintab(k-1,1));
end end end
%%%%%%%%%%%%%RIght hand side%%%%%%%%%%%%%%%%%%%%%
pCr=1; %%now start the same for the right hand side of the histogram
MM=max(find(Ensemble.peak_L{a,k1}==Ensemble.maxpeak{a,k1}));
MaP=Ensemble.maxpeak{a,k1};
k=min(MM+2,length(ZZx{a,k1})-1);
while pCr>=ccut & k<length(ZZx{a,k1})
k=min(MM+2,length(ZZx{a,k1})-1);
while k<length(ZZx{a,k1}) & Diff(k)<=0
k=k+1;
end
Minbreakright=k;
while k<=length(ZZx{a,k1})-1 & Ensemble_peak{a,k1}(k)==0
k=k+1;
end
eins=k;
k=Minbreakright;
while k<length(ZZx{a,k1})-1 & Diff(k)>=0 %%find first Max after Min to the left
k=k+1;
end if k==length(ZZx{a,k1})-1
zwei=k+1;
else if ZZy{a,k1}(k+1)>ZZy{a,k1}(k)
zwei=k+1;
else
zwei=k;
end end
if eins==length(ZZx{a,k1}) & zwei==length(ZZx{a,k1})
Nextsignifpeakright=length(ZZx{a,k1})-1;
Peaknor=Ensemble.peak_no{a,k1};
pCr=0;
Minbreakright=length(ZZx{a,k1})-1;
else
Nextsignifpeakright=min(eins,zwei);
Peaknor=-1;
for k=1:Ensemble.peak_no{a,k1}
c=find(Ensemble.peak_L{a,k1}==k);
for kk=1:length(c) if c(kk)==Nextsignifpeakright
Peaknor=k;
end end end
EE=zeros(1,length(ll));
zz=0;
if Peaknor>MaP & Peaknor~=-1 for k=MaP:Peaknor-1
EE=EE+Ensemble.signal{a,k1,k}.*L{a,k1,k};
zz=zz+L{a,k1,k};
end else
EE=Ensemble.signal2{a,k1,Nextsignifpeakright}./length(Ensemble.bars{a,Nextsignifpeakright,k1});
end if Peaknor>MaP & Peaknor~=-1
EE=EE./zz;
Sigmax=EE;
Sigright=Ensemble.signal{a,k1,Peaknor};
LUmax=median(Sigmax);
LUr=median(Sigright);
CORRCOEFr=corrcoef(Sigmax-LUmax,Sigright-LUr);
pCr=CORRCOEFr(2,1);
MM=max(find(Ensemble.peak_L{a,k1}==Peaknor));
MaP=Peaknor;
elseif Peaknor==-1
EEmax=zeros(1,length(ll));
z=0;
for k=Maxhist{a,k1}:Minbreakright
EEmax=EEmax+Ensemble.signal2{a,k1,k};
z=z+1;
end
EEmax=EEmax./z;
Sigmax=EEmax;
Sigright=EE;
LUmax=median(Sigmax);
LUr=median(Sigright);
CORRCOEFl=corrcoef(Sigmax-LUmax,Sigright-LUr);
pCr=CORRCOEFl(2,1); %%check if corresponding cluster peak is highly correlated to maximum cluster peak
c=find(ZZy{a,k1}>0);
k=Nextsignifpeakright+1;
while k<length(ZZx{a,k1}) & ZZy{a,k1}(k)==0
k=k+1;
end
MM=k;
% MaP=Nextsignifpeakright; else
pCr=0;
Peaknor=Peaknor+1;
end end end if Peaknor==Ensemble.peak_no{a,k1} & Minbreakright>=max(find(Ensemble.peak_L{a,k1}==Ensemble.peak_no{a,k1}))
cr=Minbreakright;
elseif Peaknor==1
cr=max(find(Ensemble.peak_L{a,k1}==Peaknor));
k=length(mintab(:,1));
while cr<=mintab(k,1)
k=k-1;
end
cr=min(length(mintab(:,1)),mintab(k+1,1));
elseif Peaknor==-1
cr=Minbreakright;
else
Peaknor=Peaknor-1;
% if Peaknor==Ensemble.maxpeak{a,k1}
cr=max(find(Ensemble.peak_L{a,k1}==Peaknor));
k=length(mintab(:,1));
if cr<mintab(end,1)
clear rightpeaks leftpeaks clusterbars CORRCOEFl CORRCOEFr Ensemble.signal %
%%------------------------------------------------------------------------
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%Alternative clusterbar selection by only taking those in general with a
%%%CC>0.8%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for k1=1:num
for a=1:length(data.FF{2,1})
a
ll=(max(a-rw./2+1,1):min(length(data.gpmean{2,1}),a+rw./2));
%%now the second alternative version: take signals from each
%%ind. mito:
for k=1:length(data.FF) ifisempty(data.FF{k,k1})==0
signal3{k,k1}=data.gpmean{k+1,k1}(ll);
end end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%now do the CC%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% cc1(1)=min(find(Ensemble.peak_L{a,k1}==Ensemble.maxpeak{a,k1})); % cc1(2)=max(find(Ensemble.peak_L{a,k1}==Ensemble.maxpeak{a,k1}));
% clustermaxmin=min(find(Ensemble.peak_L{a,k1}==Ensemble.maxpeak{a,k1})); % clustermaxmax=max(find(Ensemble.peak_L{a,k1}==Ensemble.maxpeak{a,k1})); % clustermax=clustermaxmin:1:clustermaxmax;
clustermax=cmmi{k1}(a):cmma{k1}(a);
% clear clustermaxmitos
clustermaxmitos{a,k1}=Ensemble.bars{a,clustermax(1),k1};
iflength(clustermax)>=2 for kk=2:length(clustermax)
clustermaxmitos{a,k1}=[clustermaxmitos{a,k1},Ensemble.bars{a,clustermax(kk),k1}];
end end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%this using the second alternative version%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%5 for k=1:length(data.FF) ifisempty(data.FF{k,k1})==0 ifisempty(find(clustermaxmitos{a,k1}==k))==1 clear Sigmax
Sigmax=Ensemble.signal{a,k1,Ensemble.maxpeak{a,k1}};
LUmax=median(Sigmax);
Sigtot=signal3{k,k1};
LUtot=median(Sigtot);
CORRCOEFtot=... corrcoef(Sigmax-LUmax,...
Sigtot-LUtot);
pCtot=CORRCOEFtot(2,1);
if pCtot>=cccut_alt2
totpeaks{a,k1}(k)=1;
else
totpeaks{a,k1}(k)=0;
end else
totpeaks{a,k1}(k)=1;
end end end clear cctot
cctot=find(totpeaks{a,k1}==1);
%--------------------------------------------------------------------------
%%%%determine max and min frequencies within the range of 1STD
%%%%forall mitos that are >=cccut_alt2 correlated to
%%%%the max peak
% clear dd % for i=1:length(cctot) % dd(i)=data.FF{cctot(i),k1}(a); % end % Ensemble.maxmean{k1}(a)=mean(dd(); % stdd=std(dd(); % % maxx{a,k1}=max(dd(); % % minn{a,k1}=min(dd(); % Ensemble.maxfrequ{k1}(a)=max(ZZx{a,k1}(2),Ensemble.maxmean{k1}(a)-stdd); % Ensemble.minfrequ{k1}(a)=min(ZZx{a,k1}(end-2),Ensemble.maxmean{k1}(a)+stdd);
dplus{a,k1}=cctot;
for i=1:length(cctot) if FFinv{a,k1}(cctot(i))==min(ZZx{a,k1}(:)) | FFinv{a,k1}(cctot(i))==max(ZZx{a,k1}(:))
cctot(i)=0;
end end
dminus{a,k1}=cctot(find(cctot>0));
ifisempty(dminus{a,k1})==1 % ctot{a,1,k1}=0; % clustertot=0;
Ensemble.clustermitos_alternative2{a,k1}=0;
else
Ensemble.clustermitos_alternative2{a,k1}=dminus{a,k1};
end
end end clear clustermitos1
clustermitos1=Ensemble.clustermitos_alternative2;
clustermitos1=dplus;
clear CORRCOEFl CORRCOEFr CORRCOEFtot pCr pCl pCtot signal2 signal3
% clear rightpeaks leftpeaks totpeaks version=3;
clear mitosinclusterr
for k1=1:num
for a=1:length(data.FF{2,1})
a
% mito_no{version,k1}(a)=length(find(Bwcluster{a,k1}==1))./length(find(data.bww1{k1}));
mitosinclusterr{version,k1}(a)=length(clustermitos1{a,k1})./length(data.FF);
clear BBvec
BBvec=max(lCUT,ZZZx{a,k1,kkk1}(1)):seg:min(ZZZx{a,k1,kkk1}(end),uCUT); %define frequ vector being the basis forall myocytes
Clusterhisty{a,k1,kkk1}=zeros(1,length(lCUT:seg:uCUT)); %predefine histogram y-scale
for b=1:length(BBvec) if cmmaa{k1,kkk1}(a)>length(BBvec) %%check that upper frequ in single histogram is within the limits of frequ vector
cmmaa{k1,kkk1}(a)=length(BBvec);
end if cmmii{k1,kkk1}(a)>cmmaa{k1,kkk1}(a) %%check that upper frequ
cmmii{k1,kkk1}(a)=cmmaa{k1,kkk1}(a);
end if BBvec(b)>=ZZZx{a,k1,kkk1}(1) & BBvec(b)<=ZZZx{a,k1,kkk1}(end) if BBvec(b)<ZZZx{a,k1,kkk1}(cmmii{k1,kkk1}(a)) | BBvec(b)>ZZZx{a,k1,kkk1}(cmmaa{k1,kkk1}(a))
Clusterhisty{a,k1,kkk1}(round((BBvec(b)-lCUT)*(1/seg))+1)=0;
else
Clusterhisty{a,k1,kkk1}(round((BBvec(b)-lCUT)*(1/seg))+1)=...
ZZZy{a,k1,kkk1}(round((BBvec(b)-ZZZx{a,k1,kkk1}(1))*(1/seg))+1);
end end end
Nsingle{k1,kkk1}=Nsingle{k1,kkk1}+mitosincluster{kkk1}{3,k1}(a).*Ntott{kkk1,k1};
Meanclusterhisty{k1,kkk1}=Meanclusterhisty{k1,kkk1}+...
Clusterhisty{a,k1,kkk1}./(mitosincluster{kkk1}{3,k1}(a).*Ntott{kkk1,k1});
z1{kkk1,k1}=z1{kkk1,k1}+1;
end end
Meanclusterhisty{k1,kkk1}=Meanclusterhisty{k1,kkk1}./z1{kkk1,k1};
end end
Meanmeanclusterhisty=zeros(1,length(lCUT:seg:uCUT));
Meanmeanclusterhistx=lCUT:seg:uCUT;
Meanmeanclusterhistx=Meanmeanclusterhistx';
Ntot=0;
A(:,1)=Meanmeanclusterhistx;
for kkk1=1:L
XX=Meanclusterhisty{1,kkk1} %%./max(Meanclusterhisty{1,kkk1}(:))
Meanmeanclusterhisty=Meanmeanclusterhisty+XX;
Ntot=Ntot+Nsingle{1,kkk1};
A(:,kkk1+1)=XX;
end
Meanmeanclusterhisty=Meanmeanclusterhisty';
figure; plot(Meanmeanclusterhisty) %
YY=Meanmeanclusterhisty;
XXX=Meanmeanclusterhistx;
YY = YY./max(YY);
ff=fittype('gauss1');
fopts=fitoptions('gauss1');
c = find(YY~=0);
xx = XXX(c);
yy = YY(c);
[cfun,gof,output]...
=fit(xx,yy,ff);
fit(Meanmeanclusterhistx,Meanmeanclusterhisty,ff);
figure;
hold on
plot(XXX,YY)
prompt={['Save files as Overall_frequs']};
defans={pathname};
fields = {'name'};
options.Resize='on';
options.WindowStyle='normal';
info = inputdlg(prompt,'File',1, defans,options);
if ~isempty(info) %see if user hit cancel
info = cell2struct(info,fields);
ii=info.name;
save([ii,'Overall_frequs_WH' ],...
'ZZZx','ZZZy','cmmii','cmmaa','Meanclusterhisty',...
'cfun','gof','output','limits',...
'Meanmeanclusterhisty','Meanmeanclusterhistx','file','bb','XXX','YY','GAUSS1','bbend','Nsingle','Ntot');
end
prompt={['Save files as Overall_frequs_origin']};
defans={ii};
fields = {'name'};
options.Resize='on';
options.WindowStyle='normal';
info = inputdlg(prompt,'File',1, defans,options);
if ~isempty(info) %see if user hit cancel
info = cell2struct(info,fields);
ii=info.name;
save([ii,'Overall_frequs_WH_origin' ],...
'A','YY','GAUSS1');
end
Angehängt findest du ausserdem die drei Files, die zum Laufen des Codes eingefügt werden soll. die geplottete und gefittete Variable ist YY.
Ich habe keine Ahnung, ob du mit diesen Dingen nun eher sagen kannst, was ich ändern muss... wenn du noch irgendwelche Informationen brauchst, lass es mich wissen! Ich danke dir sehr für deine Hilfe!
das ist sehr viel Code für ein recht kleines Problem. Wenn du nur die an deinem ursprünglich geposteten Code-Abschnitt beteiligten Variablen in eine .mat-Datei abspeicherst, würde das vollkommen reichen.
Hallo Harald,
ich verstehe, tut mir Leid Ich glaube, im Anhang sollten jetzt alle drin sein... ich habe die restlichen rausgelöscht um es übersichtlicher zu machen.
Hilft das?
LG Larissa
Verfasst am: 13.04.2017, 16:06
Titel: Programm läuft jetzt durch, fit ist aber suboptimal. Lösun
Hallo lieber Harald,
vielen Dank für deine Hilfe! Ich habe aus deinem vorherigen Kommentar nicht ersehen können, was du meintest. Danke für die idiotensichere Vorgabe eines Ersatzabschnitts der Code läuft jetzt durch, aber wie in einigen anderen Fällen ist der Fit nicht direkt optimal. Hast du Tipps, wie ich das optimieren kann?
Im Anhang der Fit zu den Daten inkl. Plot der Daten.
VLG
Larissa
Hallo Harald,
wie würde ich diese Startwerte ändern müssen im Code? Ich bin blutiger Anfänger und zudem Mediziner, kein Mathematiker, sodass ich mit statistischem Fitting nicht soooo vertraut bin. Der code ist von meinem Betreuer, sodass es durchaus sein kann, dass ein anderer Fit passender wäre... hättest du einen vorschlag?
Zum Vergleich m Anhang ein guter Fit, obwohl auch dort mehrere Peaks vorhanden sind.
wie du Schranken und Startwerte anpasst, findest du in der Doku von
fit
. Konkretere Vorschläge für die Startwerte habe ich nicht.
Was soll denn an dem angehängten Fit gut sein? Einzelne Peaks sind deutlich höher, dazwischen sind etliche Nullbereiche. Das hat doch nichts mit einer Gauß-Kurve zu tun?
Für mich sieht das nach einem Spektrum aus. In dem Fall wäre die Frage, warum man da überhaupt eine Kurve durchlegen will. Da wären doch nur die Peaks interessant?
Vielleicht mal mit dem Betreuer darüber reden...
Grüße,
Harald
Einstellungen und Berechtigungen
Du kannst Beiträge in dieses Forum schreiben. Du kannst auf Beiträge in diesem Forum antworten. Du kannst deine Beiträge in diesem Forum nicht bearbeiten. Du kannst deine Beiträge in diesem Forum nicht löschen. Du kannst an Umfragen in diesem Forum nicht mitmachen. Du kannst Dateien in diesem Forum posten Du kannst Dateien in diesem Forum herunterladen
MATLAB, Simulink, Stateflow, Handle Graphics, Real-Time Workshop, SimBiology, SimHydraulics, SimEvents, and xPC TargetBox are registered trademarks and The MathWorks, the L-shaped membrane logo, and Embedded MATLAB are trademarks of The MathWorks, Inc.
);
Ich glaube, im Anhang sollten jetzt alle drin sein... ich habe die restlichen rausgelöscht um es übersichtlicher zu machen.
RSS