vof_strengthofevidence_withplot.m

function vof_strengthofevidence_withplot(feFileToLoad, fibername, roisDir, roiNames, roiOperations, saveDir, nboots, nmontecarlo)

% Compute the Strength of Connection Evidence (S) for VOF. Creating plots
% describing the distribution of RMSE of model with and without VOF.
% Generating VOF fiber file (.pdb format).
%
% INPUT:
% feFileToLoad: a full path to .mat file containing fe structure (generated by LIFE code). Use strings for loading multiple connectome data.
% fibername: name for fiber groups. Used for file name of fiber file.
% roisDir: a directory containing waypoint ROI mat files
% roiNames: names for waypoint ROIs (string)
% roiOperations: default is {'and', 'and'} (segmenting fibers passing through both waypoint ROIs)
% saveDir:  a full path to a directory in which plot image files will be saved.
% nboots: number of bootstrap. Default is 100000.
% nmontecarlo: number of the repetition of montecarlo simulation. Default is 5.%
%
% EXAMPLE:
% feFileToLoad{1} = '/home/vof/data/S1/life/LH_1stconnectome_fe.mat';
% feFileToLoad{2} = '/home/vof/data/S1/life/LH_2ndconnectome_fe.mat';
% feFileToLoad{3} = '/home/vof/data/S1/life/LH_3rdconnectome_fe.mat';
% fibername = 'LVOF';
% roisDir = '/home/vof/data/S1/ROIs/Waypoint/'
% roiNames = {'LH_VOFPlane1.mat','LH_VOFPlane2.mat'};
% roiOperations = {'and','and'};
% saveDir = '/home/vof/data/S1/life/Picture/';
% vof_strengthofevidence_withplot(feFileToLoad, fibername, roisDir, roiNames, roiOperations, saveDir);
% 
% Copyright (C) Hiromasa Takemura and Franco Pestilli, 2013, Stanford VISTA team

% Argument checking
if notDefined('roiOperations')
    roiOperations = {'and','and'};
end
if notDefined('nboots')
    nboots = 100000;
end
if notDefined('nmontecarlo')
    nmontecarlo = 5;
end
sizecon = size(feFileToLoad);

for connum = 1:sizecon(2)
    
    % Define output file names
    fibernamenum = [fibername num2str(connum)];
    savefibersfile = [fibernamenum, '.pdb'] ;
    savematfile = [fibername ,'_stats.mat'];
    
    % Load rois
    for iroi = 1:length(roiNames)
        rois{iroi} = fullfile(roisDir,roiNames{iroi});
    end
    
    % Load fe structure
    disp('loading the LiFE structure...')
    if ischar(feFileToLoad{connum})
        fprintf('Loading %s ...\n',feFileToLoad{connum})
        load(feFileToLoad{connum});
    else
        fe  =feFileToLoad{connum};
        clear feFileToLoad;
    end
    
    % Extract the fiber group from the FE structure
    fg = feGet(fe,'fibers acpc');
    
    % Segment VOF from connectome
    fprintf('Segmenting the Tract from Connectome ...\n')
    [fgsegment, keepFascicles] = feSegmentFascicleFromConnectome(fg, rois, roiOperations, 'prob connectome');
    
    % Removing outliers from VOF fiber group
        fprintf('Excluding outliers from Fasciculus ...\n')
    [fgsegment2 keepFascicles2] = vof_gradient_removeoutlier(fe, keepFascicles,2,rois{1},rois{2});
    [fgsegment3, keepFascicle3] = mbaComputeFibersOutliers(fgsegment2,3,3,25);

    % Save VOF fibers in pdb file
    fgWrite(fgsegment3, savefibersfile);
    
    keepFascicles4 = keepFascicles2;
    keepFascicle2matrix = find(keepFascicles2);
    
    for kj = 1:length(keepFascicle3)
        if keepFascicle3(kj)==0,
            keepFascicles4(keepFascicle2matrix(kj)) = 0;
        else
        end
    end
    
    fprintf('Testing the Tract ...\n')
    [se{connum}] = feVirtualLesion(fe, keepFascicles4);
    
end
% save file
save(savematfile, 'se');

end