regnet corr xval qualify

src/MATLAB/clusterExamples/getMetaOCT.m

+function tblMeta = getMetaOCT(ROOT,eye,scanType)
+
+tblMeta = table();
+flist = dir(fullfile(ROOT,'**/*.LEVER'));
+mx = regexp({flist.name},'P(\d+)_(.+?)_(\d+-\d+-\d+)_(\d+-\d+-\d+)_(\w\w)_(.+?)_.*.LEVER','tokens');
+tblSkip = table();
+for ff = 1:length(flist)
+    nt = table();
+
+    nt.filename = {flist(ff).name};
+    nt.folder = {flist(ff).folder};
+    if isempty(mx{ff})
+        tblSkip = [tblSkip;nt];
+        continue;
+    end
+    nt.subjectID = str2double(mx{ff}{1}{1});
+    nt.scanType = mx{ff}{1}{2};
+    nt.dx = datetime(mx{ff}{1}{3},'InputFormat','MM-dd-yyyy');
+    nt.time = {mx{ff}{1}{4}};
+    nt.eye = mx{ff}{1}{5};
+    nt.scanID =  {mx{ff}{1}{6}};
+    
+    
+    if ~strcmp(nt.eye,eye) || ~contains(lower(scanType),lower(nt.scanType))
+        continue;
+    end
+    
+    tblMeta = [tblMeta;nt];
+end
+
+tblMeta = sortrows(tblMeta,{'subjectID','filename','dx'});
\ No newline at end of file

src/MATLAB/clusterExamples/goLR.m

+tic
+ROOT = '/g/leverjs/Schuman_OCT/OCT/combined';
+
+PID = 'P10010';
+flist = dir(fullfile(ROOT,['/' PID '*.LEVER']));
+mx = regexp({flist.name},'P(?<PID>\d+)_(?<scanType>.+)_(?<date>\d+-\d+-\d+)_(?<time>\d+-\d+-\d+)_(?<eye>\w\w)_(?<scanID>.+?)_.*.LEVER','names');
+tblMeta = struct2table(vertcat(mx{:}));
+tblMeta.PID = str2double(tblMeta.PID);
+tblMeta.folder = {flist.folder}';
+tblMeta.filename = {flist.name}';
+
+tblMeta = tblMeta(contains(tblMeta.scanType,'Optic'),:);
+% tblMeta = tblMeta(contains(tblMeta.scanType,'Macula'),:);
+% tblMeta = tblMeta(contains(tblMeta.eye,'OD'),:);
+
+tblQualify = readtable('../ncdVSvfmd/nyu onh matched.xlsx');
+mx = regexp(tblQualify.file_name,'P(\d+)_(.+?)_(\d+-\d+-\d+)_(\d+-\d+-\d+)_(\w\w)_(.+?)_.*.img','tokens');
+tblQualify.scanID = cellfun(@(x) x{1}{6},mx,'UniformOutput',false);
+tblQualify.PID = cellfun(@(x) str2double(x{1}{1}),mx,'UniformOutput',false);
+PID = str2double(PID(2:end));
+idx = find(cell2mat(tblQualify.PID) == PID);
+
+for i = 1:height(tblMeta)
+    idx = find(strcmp(tblQualify.scanID,tblMeta.scanID{i}));
+    if any(strcmp(tblQualify.qualification_status(idx),'No'))
+        tblMeta.bDisqualified(i) = true;
+    end
+end
+tblMeta(tblMeta.bDisqualified,:) = [];    
+p = ljsStartParallel(96);
+
+imf = RSF.getImf(tblMeta,[[0.5,8,8] ; [5,5,5]],'');
+d = [];
+parfor i = 1:length(imf)
+    dj = [];
+    for j = 1:length(imf)
+        dj(j) = SSF.ncd_ssf_volume(imf{i},imf{j});
+    end
+    d(i,:) = dj;
+end
+
+A = Cluster.Regularize(d);
+[~,Y] = Cluster.SpectralCluster(A,3);
+
+tblMeta.Y = Y;
+
+tS = tblMeta(strcmp(tblMeta.eye,'OD'),:);
+tD = tblMeta(strcmp(tblMeta.eye,'OS'),:);
+t1 = tS;
+t2 = tD;
+
+tO = tblMeta(contains(tblMeta.scanType,'Optic'),:);
+tM = tblMeta(contains(tblMeta.scanType,'Macular'),:);
+t3 = tO;
+t4 = tM;
+
+clf;hold on
+plot3(t1.Y(:,1),t1.Y(:,2),t1.Y(:,3),'r*')
+plot3(t2.Y(:,1),t2.Y(:,2),t2.Y(:,3),'og')
+legend({'OS','OD'})
+xlabel('NCD1')
+ylabel('NCD2')
+zlabel('NCD3')
+
+plot3(t3.Y(:,1),t3.Y(:,2),t3.Y(:,3),'mx')
+plot3(t4.Y(:,1),t4.Y(:,2),t4.Y(:,3),'cs')
+legend({'ONH','macula'})
+toc
\ No newline at end of file

src/MATLAB/ncdVSvfmd/getStats.m

@@ -15,12 +15,21 @@ for i = 1:length(dxx)
             continue
         end
        
+        if pData.dx(m) > pData.dx(n)
+            4;
+            tmp = m;
+            m = n;
+            n = tmp;
+        end
        if m~=n && pData.dx(m) == pData.dx(n) % same date pairs
             continue
         end
 %         
         
         dmd = pData.md(n) - pData.md(m);
+        if dmd>0
+            continue
+        end
         nt = table();
         nt.subjectID = pData.subjectID(m);
         idxm = find(dxx{i}(:,2)==m & dxx{i}(:,3) == m);
@@ -28,18 +37,18 @@ for i = 1:length(dxx)
         nt.X = [dxx{i}(idxm,1),dxx{i}(idxn,1),dxx{i}(j,1)];
         nt.velocity = (dxx{i}(j,1));        
         nt.md = min(pData.md(n),pData.md(m));
-        nt.dmd = abs(dmd);
+        nt.dmd = dmd;
         nt.mn = [m,n];                
         nt.ddx = pData.dx(n) - pData.dx(m);
         % find sex
         idx = find(tblMD.subjectID == pData.subjectID(m),1,'first');
         nt.sex = tblMD.gender_full(idx);
-        if nt.dmd > 4
-             continue
-        end
-        if years(nt.ddx)>1
+        if abs(nt.dmd)>4
             continue
         end
+%         if years(nt.ddx)>1
+%             continue
+%         end
 %         if nt.velocity > 0.9 && nt.velocity < 0.975
 %             continue
 %         end

src/MATLAB/ncdVSvfmd/goQualify.m

 load('onh_kernel_corr_pt5_5_15.mat')
+% load('kcorr_VFMD_plate_blob.mat')
+% load('kernel_corr_raw.mat')
 tblQualify = readtable('nyu onh matched.xlsx');
 mx = regexp(tblQualify.file_name,'P(\d+)_(.+?)_(\d+-\d+-\d+)_(\d+-\d+-\d+)_(\w\w)_(.+?)_.*.img','tokens');
 tblQualify.scanID = cellfun(@(x) x{1}{6},mx,'UniformOutput',false);
-
+tblQualify.patientID = cellfun(@(x) str2double(x{1}{1}),mx,'UniformOutput',false);
 for i = 1:height(pData)
 
     idxScan = find(strcmp(tblQualify.scanID,pData.scanID{i})); 
@@ -14,3 +16,5 @@ for i = 1:height(pData)
 end
 
 
+tblVelo = getStats(dxx,pData,tblVFMD);
+[err,mdlX] = regnetCorr(tblVelo);

src/MATLAB/ncdVSvfmd/regnetCorr.m

-function [rmse1,rmse2] = regnetCorr(tblVelo)
+function [err,mdlX] = regnetCorr(tblVelo)
+
+% tblVelo.velocity is NCD 
+% tblVelo.md is vfmd for first of two image pairs
+% tblVelo.dmd is delta VFMD between two image pairs
+
 % tblVelo = getStats(dxx,pData,[],tblVFMD);
 % X = [tblVelo.velocity,tblVelo.md];
 t2 = tblVelo(tblVelo.dmd~=0,:);
@@ -7,22 +12,23 @@ X1 = [t2.velocity];
 X2 = [t2.velocity,t2.md];
 T = t2.dmd;
 
-% mdl1 = fitrnet(X1,T,'Activations','relu','LayerSizes',[44]);
-% cv1 = mdl1.crossval; rmse1 = cv1.kfoldLoss^0.5;
-% mdl2 = fitrnet(X2,T,'Activations','sigmoid','LayerSizes',[10]);
-% cv2 = mdl2.crossval; rmse2 = cv2.kfoldLoss^0.5;
-% 4;
-% mdl = fitrnet(X,T,"OptimizeHyperparameters","auto", "HyperparameterOptimizationOptions",struct("AcquisitionFunctionName","expected-improvement-plus"))
-
+% leave one out x-validation
 err = [];
-for i = 1:length(X1)
-    trainX2 = X1;
-    trainX2(i) = [];
-    trainY2 = T;
-    trainY2(i) = [];
-    mdl = fitrnet(trainX2,trainY2,'Activations','sigmoid','LayerSizes',[10]);
-    pred = mdl.predict(X2(i));
+mdlX = {};
+xTarget = X1;
+for i = 1:length(xTarget)
+    trainX = xTarget;
+    trainX(i,:) = [];
+    trainY = T;
+    trainY(i) = [];
+%     mdl = fitrnet(trainX2,trainY2,'Activations','sigmoid','standardize',false,'LayerSizes',[10]);
+    mdl = fitrnet(trainX,trainY,"OptimizeHyperparameters","auto", "HyperparameterOptimizationOptions",...
+        struct("AcquisitionFunctionName","expected-improvement-plus",'verbose',0,'showplots',false));
+    mdlX{i} = mdl;
+    pred = mdl.predict(xTarget(i,:));
     err(i) = abs(T(i) - pred);
+    [i,mean(err)]
 end
+mean(err)
 4;