Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

add outlier functions to bypass matlab's isoutlier; adapt to r2017b #263

Merged
merged 9 commits into from
May 9, 2023
Merged

Conversation

cyang31
Copy link
Contributor

@cyang31 cyang31 commented Apr 25, 2023

Background

Matlab's isoutlier function has new option percentile that wasn't introduced in r2017b but in the latest versions. Because of this, the compiled version of Matlab always throws errors for scripts written in r2022 which use commands such as isoutlier(a, 'percentile', [5,95]).

Main changes

  • isoutlier and filloutlier are removed from the tICA feature generation step
  • two custom outlier detection functions are involved,
    • findFivePercentageOutliers to find outliers outside 5% to 95% range
    • findScaledMADOutliers to find outliers three scaled MAD away from the median
  • the final outlier is determined by the intersection (AND) between the outliers obtained from the above criteria

Evaluation

In r2022b (raw version),

A = [57 59 60 100 59 58 57 58 300 61 62 60 62 58 57];
res=isoutlier(A,'percentile',[5,95])

res =

  1×15 logical array

   0   0   0   0   0   0   0   0   1   0   0   0   0   0   0

in r2022b & r2017b (current commit),

A = [57 59 60 100 59 58 57 58 300 61 62 60 62 58 57];
res=findOutliers(A)

res =

  1×15 logical array

   0   0   0   0   0   0   0   0   1   0   0   0   0   0   0

@cyang31 cyang31 requested review from coalsont and glasserm April 25, 2023 03:38
@coalsont
Copy link
Member

To be clear, the "inconsistency" is just that the percentile option was added after r2017b, right? I wouldn't call that "inconsistent", it is just a new feature that we didn't realize was missing from the older version of the MCR we use.

@cyang31
Copy link
Contributor Author

cyang31 commented Apr 25, 2023

Fair enough.

@mharms
Copy link
Contributor

mharms commented Apr 25, 2023

Should we consider updating the MCR version that the compiled matlab uses? Or do we not want to go there...?

@cyang31
Copy link
Contributor Author

cyang31 commented Apr 25, 2023

Should we consider updating the MCR version that the compiled matlab uses? Or do we not want to go there...?

I think the reason why we keep using the r2017b is because the non-python version of FIX is no longer being updated and it is using r2017b runtime. Perhaps this raises another question if we should switch to the pyFIX at some point in the future...

@mharms
Copy link
Contributor

mharms commented Apr 25, 2023

Yes, the container and related environment variables would need to support two different MCRs. But that was in fact the situation at one point previously, so its doable.

@mharms
Copy link
Contributor

mharms commented Apr 25, 2023

At some point, given the amount of Matlab we have other than FIX, we are going to want to update the Matlab MCR version, although maybe it does make sense to wait to do that (if we can wait) until we switch over to using pyFIX.

@coalsont
Copy link
Member

This looks reasonable to me. Any other opinions?

@glasserm
Copy link
Contributor

I am okay with committing if this is tested and Tim approves it.

@cyang31
Copy link
Contributor Author

cyang31 commented Apr 28, 2023

I have tested it and the new feature csv is generated here without error: /media/myelin/brainmappers/Connectome_Project/YA_HCP_Final/S1200_MSMAll3T535T/MNINonLinear/Results/rfMRI_REST/tICA_d84/features.csv

@coalsont
Copy link
Member

I'll need to recompile it against r2017b, preferably before merging.

@cyang31
Copy link
Contributor Author

cyang31 commented Apr 28, 2023

If you want me to test the compiled version, please let me know when it's ready.

@coalsont
Copy link
Member

I have pushed a commit with that function recompiled, try that in compiled mode.

@cyang31
Copy link
Contributor Author

cyang31 commented Apr 28, 2023

I guess I just need to git pull to get the new script?

@coalsont
Copy link
Member

Yes.

@cyang31
Copy link
Contributor Author

cyang31 commented Apr 28, 2023

Hi Tim, I get an error

LD_LIBRARY_PATH is .:/export/matlab/MCR/R2017b/v93/runtime/glnxa64:/export/matlab/MCR/R2017b/v93/bin/glnxa64:/export/matlab/MCR/R2017b/v93/sys/os/glnxa64:/export/matlab/MCR/R2017b/v93/sys/opengl/lib/glnxa64 /media/myelin/alex/pipelines/pipeline1/HCPpipelines/tICA/scripts/Compiled_ComputeTICAFeatures/ComputeTICAFeatures: error while loading shared libraries: libmwlaunchermain.so: cannot open shared object file: No such file or directory

The error log is saved here /media/myelin/alex/tICA/tICAClassify/tICAClassifier/prepare_data/error.log

@coalsont
Copy link
Member

I would guess you have the MCR set incorrectly in your SetUp... script, /export isn't where we put stuff on our lab machines. I don't know if we have MCR v93 installed on them.

@cyang31
Copy link
Contributor Author

cyang31 commented Apr 28, 2023

Oh, I see. There's nothing under this path /export/matlab/MCR/R2017b/v93. Is it quick to install one on desktop4 such that I can test it for the compiled version?

@coalsont
Copy link
Member

Try /mnt/myelin/tim/matlabinst/MCR_install/v93.

@cyang31
Copy link
Contributor Author

cyang31 commented Apr 28, 2023

Here's another error.

Can't reload '/media/myelin/tim/matlabinst/MCR_install/v93/bin/glnxa64/libmwcoder_types.so'

Error in aryule (line 71)


Error in single_subject_tcs_features (line 21)


Error in ComputeTICAFeatures (line 623)

MATLAB:dispatcher:loadLibrary
Fri Apr 28 18:01:49 CDT 2023:ComputeTICAFeatures.sh: While running '/media/myelin/alex/pipelines/pipeline1/HCPpipelines/tICA/scripts/ComputeTICAFeatures.sh --study-folder=/media/myelin/brainmappers/Connectome_Project/YA_HCP_Final --out-group-name=S1200_MSMAll7T87T --subject-list=102311@105923@109123@111312@111514@114823@115017@125525@126426@128935@130518@132118@134829@137128@144226@145834@146129@146432@146735@146937@148133@157336@158035@159239@164636@167440@169343@169444@169747@175237@176542@177746@178142@178647@180533@182436@185442@186949@191033@192439@195041@196144@197348@200311@200614@204521@205220@214019@214524@221319@246133@249947@263436@320826@346137@380036@381038@393247@406836@429040@436845@541943@550439@562345@572045@585256@601127@627549@638049@654552@690152@706040@732243@751550@771354@782561@783462@818859@859671@861456@872764@899885@901442@910241@942658@966975@995174 --fmri-list=tfMRI_RETCCW_7T_AP@tfMRI_RETCW_7T_PA@tfMRI_RETEXP_7T_AP@tfMRI_RETCON_7T_PA@tfMRI_RETBAR1_7T_AP@tfMRI_RETBAR2_7T_PA --fmri-output-name=tfMRI_RET_7T --ica-dim=75 --proc-string=_hp0_clean_rclean --tica-proc-string=tfMRI_RET_7T_d75_WF5_S1200_MSMAll7T87T_WR_tICA --fmri-resolution=1.60 --surf-reg-name=MSMAll --low-res=32 --melodic-high-pass=0 --mrfix-concat-name=tfMRI_7T_RETCCW_AP_RETCW_PA_RETEXP_AP_RETCON_PA_RETBAR1_AP_RETBAR2_PA --reclean-mode=YES --matlab-run-mode=0':
Fri Apr 28 18:01:49 CDT 2023:ComputeTICAFeatures.sh: ERROR: '/media/myelin/alex/pipelines/pipeline1/HCPpipelines/tICA/scripts/Compiled_ComputeTICAFeatures/run_ComputeTICAFeatures.sh' command failed with return code: 255

===> ERROR: Command returned with nonzero exit code
---------------------------------------------------
         script: ComputeTICAFeatures.sh
stopped at line: 140
           call: "${matlab_cmd[@]}"
  expanded call: /media/myelin/alex/pipelines/pipeline1/HCPpipelines/tICA/scripts/Compiled_ComputeTICAFeatures/run_ComputeTICAFeatures.sh /mnt/myelin/tim/matlabinst/MCR_install/v93 /media/myelin/brainmappers/Connectome_Project/YA_HCP_Final S1200_MSMAll7T87T /media/myelin/brainmappers/Connectome_Project/YA_HCP_Final/S1200_MSMAll7T87T/MNINonLinear/Results/tfMRI_RET_7T/tICA_d75/SubjectList.txt /media/myelin/brainmappers/Connectome_Project/YA_HCP_Final/S1200_MSMAll7T87T/MNINonLinear/Results/tfMRI_RET_7T/tICA_d75/fMRIList.txt tfMRI_RET_7T 75 _hp0_clean_rclean tfMRI_RET_7T_d75_WF5_S1200_MSMAll7T87T_WR_tICA 1.60 _MSMAll 32 YES 0 tfMRI_7T_RETCCW_AP_RETCW_PA_RETEXP_AP_RETCON_PA_RETBAR1_AP_RETBAR2_PA YES /media/myelin/alex/pipelines/pipeline0/HCPpipelines/global/config/tICA /media/myelin/alex/pipelines/pipeline1/HCPpipelines/tICA/scripts/feature_helpers /media/myelin/alex/pipelines/pipeline0/HCPpipelines/global/templates/tICA/Q1-Q6_RelatedValidation210.CorticalAreas_dil_Final_Final_Areas_Group_Colors.32k_fs_LR.dlabel.nii /media/myelin/alex/pipelines/pipeline0/HCPpipelines/global/templates/tICA/rfMRI_REST_Atlas_MSMAll_2_d41_WRN_DeDrift_hp2000_clean_tclean_nobias_vn_BC_CorticalAreas_dil_210V_MPM_Group_group_z_mean_TestII.txt /media/myelin/alex/pipelines/pipeline0/HCPpipelines/global/templates/tICA/Nifti_Template.1.60.nii.gz /media/myelin/alex/pipelines/pipeline0/HCPpipelines/global/templates/tICA/Vascular_Territory.1.60.nii.gz /media/myelin/alex/pipelines/pipeline0/HCPpipelines/global/templates/tICA/Vessel_Probabilities.1.60.nii.gz /media/myelin/alex/pipelines/pipeline0/HCPpipelines/global/templates/tICA/Multiband_Kspace.mat /media/myelin/alex/pipelines/pipeline0/HCPpipelines/global/templates/tICA/Partial.pvcorr_perfusion_calib_Atlas.dscalar.nii /media/myelin/alex/pipelines/pipeline0/HCPpipelines/global/templates/tICA/Partial.arrival_Atlas.dscalar.nii
       hostname: brainmappers-desktop4
      exit code: 255
---------------------------------------------------

===> Aborting execution!

@cyang31
Copy link
Contributor Author

cyang31 commented Apr 28, 2023

Hi Tim, I made a minor change to the tICA/scripts/ComputeTICAFeatures.m because it used to be too slow to run this step on a pre-computed dataset.

@coalsont
Copy link
Member

I haven't seen that "can't reload" error before. I could install the MCR locally instead of on nfs, and see if that helps.

@coalsont
Copy link
Member

Try /home/brainmappers/MCR_install/v93 on desktop4. I am recompiling now.

@cyang31
Copy link
Contributor Author

cyang31 commented Apr 28, 2023

I will try to run after I get the updated compiled script, given that the old version takes about 1hr to run to where the error occurs.

@cyang31
Copy link
Contributor Author

cyang31 commented May 3, 2023

I tested the rewritten version of aryule using the following script,

rng(0)
ss_tcs_data=rand(1000,1);
var_all=zeros(1,8);
for j_order=1:8
    [~,e,~]=aryule(double(ss_tcs_data'),j_order);
    var_all(j_order)=e;
end

var_all_myfunc=zeros(1,8);
for j_order=1:8
    [~,e]=myfunc_aryule(ss_tcs_data,j_order);
    var_all_myfunc(j_order)=e;
end

ar1=aryule(ss_tcs_data,1);

ar1_myfunc=myfunc_aryule(ss_tcs_data,1);

ar2=aryule(ss_tcs_data,2);

ar2_myfunc=myfunc_aryule(ss_tcs_data,2);

Here're the results,

>> var_all

var_all =

    0.1381    0.1117    0.1054    0.0998    0.0969    0.0936    0.0925    0.0910

>> var_all_myfunc

var_all_myfunc =

    0.1381    0.1117    0.1054    0.0998    0.0969    0.0936    0.0925    0.0910

>> ar1

ar1 =

    1.0000   -0.7532

>> ar1_myfunc

ar1_myfunc =

   -0.7532

>> ar2

ar2 =

    1.0000   -0.4240   -0.4370

>> ar2_myfunc

ar2_myfunc =

   -0.4240   -0.4370

So, I think this rewritten myfunc_aryule.m should work.

@coalsont
Copy link
Member

coalsont commented May 4, 2023

Have you tested tICA in compiled mode with it?

@cyang31
Copy link
Contributor Author

cyang31 commented May 4, 2023

Unfortunately, it seems that it have solved the aryule issue but there's another issue regarding the function pwelch in signal processing toolbox. Let me look into this to see if it can be rewritten as well.

Error using signal.internal.sigcheckfloattype
Can't reload '/home/brainmappers/matlab2017b/bin/glnxa64/libmwxcp_dwarf.so'

Error in signal.internal.sigcasttofloat (line 36)


Error in gencoswin (line 15)


Error in hamming (line 22)


Error in welchparse>segment_info (line 199)


Error in welchparse (line 32)


Error in welch (line 42)


Error in pwelch (line 158)


Error in pwelch_ps (line 4)


Error in ComputeTICAFeatures (line 623)

MATLAB:dispatcher:loadLibrary

@cyang31
Copy link
Contributor Author

cyang31 commented May 4, 2023

This is weird. I wrote this part of code based from https://github.com/Washington-University/HCPpipelines/blob/master/global/matlab/nets_spectra/nets_spectra.m#L37 and now it gives me an error.

@coalsont
Copy link
Member

coalsont commented May 4, 2023

I compiled your aryule test and ran it on my linux machine, and it worked. I then tried it on desktop5, which failed with the reload error. So, it seems there is something about the OS version or packages that is causing this reload error.

@cyang31
Copy link
Contributor Author

cyang31 commented May 4, 2023

What a relief! (especially after I found pwelch is hard to rewrite...)

Is it the raw aryule or the myfunc_aryule that you have tested on your machine?

@coalsont
Copy link
Member

coalsont commented May 4, 2023

Matlab's original aryule.

@coalsont
Copy link
Member

coalsont commented May 4, 2023

Qunex is currently using CentOS 7, I believe, which is older that what we have on desktop5, so it is possible compiled things will work in it without even replacing aryule. However, this implies that tICA in compiled mode will stay broken for modern linux until we use something newer than MCR v9.3 (and any other matlab in the pipelines that uses something more complicated than pinv may also already be broken for modern linux).

@cyang31
Copy link
Contributor Author

cyang31 commented May 4, 2023

Regarding this PR, I think we should go back to f596d68 and compile it on your machine, then merge, such that the user can use it right away.

Updating the MCR seems to be a necessary step and we should discuss with @glasserm about this given that it requires a lot of changes...

@cyang31
Copy link
Contributor Author

cyang31 commented May 8, 2023

Since the user has tested it without error, I would suggest we can merge this PR now.

@glasserm
Copy link
Contributor

glasserm commented May 9, 2023

I am okay with that if Tim approves.

@coalsont coalsont merged commit 005a096 into Washington-University:master May 9, 2023
@cyang31 cyang31 deleted the fix-tica branch June 6, 2023 07:05
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

4 participants