diff --git a/.github/workflows/comment-bot.yml b/.github/workflows/comment-bot.yml index 4451632..b44ee7b 100644 --- a/.github/workflows/comment-bot.yml +++ b/.github/workflows/comment-bot.yml @@ -29,6 +29,7 @@ jobs: post({ owner: context.repo.owner, repo: context.repo.repo, comment_id: context.payload.comment.id, content: "eyes"}) + github-token: ${{ secrets.GH_TOKEN }} - name: Tag Commit run: | git clone https://${GITHUB_TOKEN}@github.com/${GITHUB_REPOSITORY} repo @@ -48,3 +49,4 @@ jobs: post({ owner: context.repo.owner, repo: context.repo.repo, comment_id: context.payload.comment.id, content: "rocket"}) + github-token: ${{ secrets.GH_TOKEN }} diff --git a/.github/workflows/test.yml b/.github/workflows/test.yml index 7fd515b..f9a0ce5 100644 --- a/.github/workflows/test.yml +++ b/.github/workflows/test.yml @@ -57,8 +57,6 @@ jobs: - if: startsWith(matrix.python, '2') run: pytest - run: codecov - env: - CODECOV_TOKEN: ${{ secrets.CODECOV_TOKEN }} matlab: if: github.event_name != 'pull_request' || github.head_ref != 'devel' runs-on: [self-hosted, python, matlab] @@ -72,8 +70,6 @@ jobs: - run: pip install -U .[dev] - run: pytest --durations-min=1 - run: codecov - env: - CODECOV_TOKEN: ${{ secrets.CODECOV_TOKEN }} - name: Post Run setup-python run: setup-python -p3.7 -Dr if: ${{ always() }} @@ -89,37 +85,23 @@ jobs: - id: dist uses: casperdcl/deploy-pypi@v2 with: - build: true + pip: true gpg_key: ${{ secrets.GPG_KEY }} password: ${{ secrets.PYPI_TOKEN }} - upload: ${{ github.event_name == 'push' && startsWith(github.event.ref, 'refs/tags') }} - - name: Changelog - run: git log --pretty='format:%d%n- %s%n%b---' $(git tag --sort=v:refname | tail -n2 | head -n1)..HEAD > _CHANGES.md - - if: github.event_name == 'push' && startsWith(github.event.ref, 'refs/tags') - id: create_release - uses: actions/create-release@v1 + upload: ${{ github.event_name == 'push' && startsWith(github.ref, 'refs/tags') }} + - id: meta + name: Changelog + run: | + echo ::set-output name=tag::${GITHUB_REF#refs/tags/} + git log --pretty='format:%d%n- %s%n%b---' $(git tag --sort=v:refname | tail -n2 | head -n1)..HEAD > _CHANGES.md + - if: github.event_name == 'push' && startsWith(github.ref, 'refs/tags') + uses: softprops/action-gh-release@v1 env: GITHUB_TOKEN: ${{ secrets.GH_TOKEN }} with: - tag_name: ${{ github.ref }} - release_name: spm12 ${{ github.ref }} beta + name: spm12 ${{ steps.meta.outputs.tag }} beta body_path: _CHANGES.md draft: true - - if: github.event_name == 'push' && startsWith(github.event.ref, 'refs/tags') - uses: actions/upload-release-asset@v1 - env: - GITHUB_TOKEN: ${{ secrets.GH_TOKEN }} - with: - upload_url: ${{ steps.create_release.outputs.upload_url }} - asset_path: dist/${{ steps.dist.outputs.whl }} - asset_name: ${{ steps.dist.outputs.whl }} - asset_content_type: application/zip - - if: github.event_name == 'push' && startsWith(github.event.ref, 'refs/tags') - uses: actions/upload-release-asset@v1 - env: - GITHUB_TOKEN: ${{ secrets.GH_TOKEN }} - with: - upload_url: ${{ steps.create_release.outputs.upload_url }} - asset_path: dist/${{ steps.dist.outputs.whl_asc }} - asset_name: ${{ steps.dist.outputs.whl_asc }} - asset_content_type: text/plain + files: | + dist/${{ steps.dist.outputs.whl }} + dist/${{ steps.dist.outputs.whl_asc }} diff --git a/.pre-commit-config.yaml b/.pre-commit-config.yaml index 2fcbb14..8d77ba8 100644 --- a/.pre-commit-config.yaml +++ b/.pre-commit-config.yaml @@ -2,7 +2,7 @@ default_language_version: python: python3 repos: - repo: https://github.com/pre-commit/pre-commit-hooks - rev: v3.4.0 + rev: v4.0.1 hooks: - id: check-added-large-files - id: check-case-conflict @@ -26,19 +26,20 @@ repos: exclude: ^(.pre-commit-config.yaml|.github/workflows/test.yml)$ args: [-i] - repo: https://gitlab.com/pycqa/flake8 - rev: 3.8.4 + rev: 3.9.2 hooks: - id: flake8 + args: [-j8] additional_dependencies: - flake8-bugbear - flake8-comprehensions - flake8-debugger - flake8-string-format - repo: https://github.com/psf/black - rev: 19.10b0 + rev: 21.7b0 hooks: - id: black - repo: https://github.com/PyCQA/isort - rev: 5.7.0 + rev: 5.9.3 hooks: - id: isort diff --git a/setup.cfg b/setup.cfg index 89f3898..3a4310a 100644 --- a/setup.cfg +++ b/setup.cfg @@ -63,6 +63,8 @@ demo= [options.entry_points] console_scripts = spm12=spm12.cli:main +[options.packages.find] +exclude=tests [options.package_data] *=*.md, *.rst, *.m diff --git a/spm12/coreg_spm_m.m b/spm12/amypad_coreg.m similarity index 87% rename from spm12/coreg_spm_m.m rename to spm12/amypad_coreg.m index 98aabcb..ae8199d 100644 --- a/spm12/coreg_spm_m.m +++ b/spm12/amypad_coreg.m @@ -1,4 +1,4 @@ -function [M, x] = coreg_spm_m(imref, imflo, costfun, sep, tol, fwhm, params, graphics, visual) +function [M, x] = amypad_coreg(imref, imflo, costfun, sep, tol, fwhm, params, graphics, visual) if visual>0 Fgraph = spm_figure('GetWin','Graphics'); Finter = spm_figure('GetWin','Interactive'); diff --git a/spm12/amypad_coreg_modify_affine.m b/spm12/amypad_coreg_modify_affine.m new file mode 100644 index 0000000..569f74c --- /dev/null +++ b/spm12/amypad_coreg_modify_affine.m @@ -0,0 +1,7 @@ +function out = amypad_coreg_modify_affine(imflo, M) + VF = strcat(imflo,',1'); + MM = zeros(4,4); + MM(:,:) = spm_get_space(VF); + spm_get_space(VF, M\MM(:,:)); + out = 0; +end diff --git a/spm12/amypad_normw.m b/spm12/amypad_normw.m new file mode 100644 index 0000000..f7570a5 --- /dev/null +++ b/spm12/amypad_normw.m @@ -0,0 +1,10 @@ +function out = amypad_normw(def_file, flist4norm) + job.subj.def = {def_file}; + job.subj.resample = flist4norm; + job.woptions.bb = [NaN, NaN, NaN; NaN, NaN, NaN]; + job.woptions.vox = [2, 2, 2]; + job.woptions.interp = 4; + job.woptions.prefix = 'w'; + spm_run_norm(job); + out=0; +end diff --git a/spm12/resample_spm_m.m b/spm12/amypad_resample.m similarity index 85% rename from spm12/resample_spm_m.m rename to spm12/amypad_resample.m index e58bc8c..c942f85 100644 --- a/spm12/resample_spm_m.m +++ b/spm12/amypad_resample.m @@ -1,4 +1,4 @@ -function out = resample_spm_m(imref, imflo, M, f_mask, f_mean, f_interp, f_which, f_prefix) +function out = amypad_resample(imref, imflo, M, f_mask, f_mean, f_interp, f_which, f_prefix) %-Reslicing parameters rflags.mask = f_mask; rflags.mean = f_mean; diff --git a/spm12/amypad_seg.m b/spm12/amypad_seg.m new file mode 100644 index 0000000..42bf585 --- /dev/null +++ b/spm12/amypad_seg.m @@ -0,0 +1,46 @@ +function [param,invdef,fordef] = amypad_seg(f_mri, spm_path, nat_gm, nat_wm, nat_csf, store_fwd, store_inv, visual) + job.channel.vols = {strcat(f_mri,',1')}; + job.channel.biasreg = 0.001; + job.channel.biasfwhm = 60; + job.channel.write = [0, 0]; + job.tissue(1).tpm = {[spm_path, filesep, 'tpm', filesep, 'TPM.nii,1']}; + job.tissue(1).ngaus = 1; + job.tissue(1).native = [nat_gm, 0]; + job.tissue(1).warped = [0, 0]; + job.tissue(2).tpm = {[spm_path, filesep, 'tpm', filesep, 'TPM.nii,2']}; + job.tissue(2).ngaus = 1; + job.tissue(2).native = [nat_wm, 0]; + job.tissue(2).warped = [0, 0]; + job.tissue(3).tpm = {[spm_path, filesep, 'tpm', filesep, 'TPM.nii,3']}; + job.tissue(3).ngaus = 2; + job.tissue(3).native = [nat_csf, 0]; + job.tissue(3).warped = [0, 0]; + job.tissue(4).tpm = {[spm_path, filesep, 'tpm', filesep, 'TPM.nii,4']}; + job.tissue(4).ngaus = 3; + job.tissue(4).native = [0, 0]; + job.tissue(4).warped = [0, 0]; + job.tissue(5).tpm = {[spm_path, filesep, 'tpm', filesep, 'TPM.nii,5']}; + job.tissue(5).ngaus = 4; + job.tissue(5).native = [0, 0]; + job.tissue(5).warped = [0, 0]; + job.tissue(6).tpm = {[spm_path, filesep, 'tpm', filesep, 'TPM.nii,6']}; + job.tissue(6).ngaus = 2; + job.tissue(6).native = [0, 0]; + job.tissue(6).warped = [0, 0]; + job.warp.mrf = 1; + job.warp.cleanup = 1; + job.warp.reg = [0, 0.001, 0.5, 0.05, 0.2]; + job.warp.affreg = 'mni'; + job.warp.fwhm = 0; + job.warp.samp = 3; + job.warp.write = [store_fwd, store_inv]; + if visual>0 + Finter = spm_figure('GetWin','Interactive'); + end + spm_jobman('initcfg'); + segout = spm_preproc_run(job); + param = segout.param{1}; + invdef = segout.invdef{1}; + fordef = segout.fordef{1}; + %disp(segout); +end diff --git a/spm12/regseg.py b/spm12/regseg.py index ca5ce31..5d5983b 100644 --- a/spm12/regseg.py +++ b/spm12/regseg.py @@ -1,6 +1,7 @@ import errno import logging import os +import re import shutil from textwrap import dedent @@ -9,12 +10,28 @@ from miutil import create_dir, hasext from miutil.imio import nii -from .utils import ensure_spm +from .utils import ensure_spm, spm_dir __author__ = ("Pawel J. Markiewicz", "Casper O. da Costa-Luis") log = logging.getLogger(__name__) +def move_files(fin, opth): + """ + Move input file path fin to the output folder opth. + """ + fdst = os.path.join(opth, os.path.basename(fin)) + shutil.move(fin, fdst) + return fdst + + +def glob_match(pttrn, pth): + """ + glob with regular expressions + """ + return (os.path.join(pth, f) for f in os.listdir(pth) if re.match(pttrn, f)) + + def fwhm2sig(fwhm, voxsize=2.0): return fwhm / (voxsize * (8 * np.log(2)) ** 0.5) @@ -25,14 +42,13 @@ def smoothim(fim, fwhm=4, fout=""): """ imd = nii.getnii(fim, output="all") imsmo = ndi.filters.gaussian_filter( - imd["im"], fwhm2sig(fwhm, voxsize=imd["voxsize"]), mode="mirror" + imd["im"], fwhm2sig(fwhm, voxsize=imd["voxsize"]), mode="constant" ) - if not fout: - fout = "{f[0]}{f[1]}_smo{fwhm}{f[2]}".format( - f=nii.file_parts(fim), fwhm=str(fwhm).replace(".", "-") + f = nii.file_parts(fim) + fout = os.path.join( + f[0], "{}_smo{}{}".format(f[1], str(fwhm).replace(".", "-"), f[2]) ) - nii.array2nii( imsmo, imd["affine"], @@ -44,7 +60,6 @@ def smoothim(fim, fwhm=4, fout=""): ), flip=imd["flip"], ) - return {"im": imsmo, "fim": fout, "fwhm": fwhm, "affine": imd["affine"]} @@ -68,7 +83,18 @@ def coreg_spm( del_uncmpr=True, save_arr=True, save_txt=True, + modify_nii=False, ): + """Rigid body registration using SPM coreg function. + Args: + imref: reference image + imflo: floating image + fwhm_ref: FWHM of the smoothing kernel for the reference image + fwhm_flo: FWHM of the smoothing kernel for the reference image + modify_nii: modifies the affine of the NIfTI file of the floating + image according to the rigid body transformation. + """ + out = {} # output dictionary sep = sep or [4, 2] tol = tol or [ 0.0200, @@ -96,7 +122,7 @@ def coreg_spm( log.debug("output path:%s", opth) create_dir(opth) - # > decompress ref image as necessary + # decompress ref image as necessary if hasext(imref, "gz"): imrefu = nii.nii_ugzip(imref, outpath=opth) else: @@ -106,8 +132,7 @@ def coreg_spm( if fwhm_ref > 0: smodct = smoothim(imrefu, fwhm_ref) - - # > delete the previous version (non-smoothed) + # delete the previous version (non-smoothed) os.remove(imrefu) imrefu = smodct["fim"] @@ -117,7 +142,7 @@ def coreg_spm( ) ) - # > floating + # floating if hasext(imflo, "gz"): imflou = nii.nii_ugzip(imflo, outpath=opth) else: @@ -127,8 +152,13 @@ def coreg_spm( if fwhm_flo > 0: smodct = smoothim(imflou, fwhm_flo) - # > delete the previous version (non-smoothed) - os.remove(imflou) + # delete the previous version (non-smoothed) + if not modify_nii: + os.remove(imflou) + else: + # save the uncompressed and unsmoothed version + imflou_ = imflou + imflou = smodct["fim"] log.info( @@ -137,10 +167,10 @@ def coreg_spm( ) ) - # run the matlab SPM coregistration + # run the MATLAB SPM registration import matlab as ml - Mm, xm = eng.coreg_spm_m( + Mm, xm = eng.amypad_coreg( imrefu, imflou, costfun, @@ -153,6 +183,11 @@ def coreg_spm( nargout=2, ) + # modify the affine of the floating image (as usually done in SPM) + if modify_nii: + eng.amypad_coreg_modify_affine(imflou_, Mm) + out["freg"] = imflou_ + # get the affine matrix M = np.array(Mm._data.tolist()) M = M.reshape(4, 4).T @@ -167,7 +202,6 @@ def coreg_spm( create_dir(os.path.join(opth, "affine-spm")) - # --------------------------------------------------------------------------- if fname_aff == "": if pickname == "ref": faff = os.path.join( @@ -181,30 +215,26 @@ def coreg_spm( "affine-spm", "affine-flo-" + nii.file_parts(imflo)[1] + fcomment + ".npy", ) - else: - - # > add '.npy' extension if not in the affine output file name + # add '.npy' extension if not in the affine output file name if not fname_aff.endswith(".npy"): fname_aff += ".npy" faff = os.path.join(opth, "affine-spm", fname_aff) - # --------------------------------------------------------------------------- - # > safe the affine transformation + # save the affine transformation if save_arr: np.save(faff, M) if save_txt: faff = os.path.splitext(faff)[0] + ".txt" np.savetxt(faff, M) - return { - "affine": M, - "faff": faff, - "rotations": x[3:], - "translations": x[:3], - "matlab_eng": eng, - } + out["affine"] = M + out["faff"] = faff + out["rotations"] = x[3:] + out["translations"] = x[:3] + out["matlab_eng"] = eng + return out def resample_spm( @@ -245,7 +275,7 @@ def resample_spm( log.debug("output path:%s", opth) create_dir(opth) - # > decompress if necessary + # decompress if necessary if hasext(imref, "gz"): imrefu = nii.nii_ugzip(imref, outpath=opth) else: @@ -253,7 +283,7 @@ def resample_spm( imrefu = os.path.join(opth, fnm) shutil.copyfile(imref, imrefu) - # > floating + # floating if hasext(imflo, "gz"): imflou = nii.nii_ugzip(imflo, outpath=opth) else: @@ -280,7 +310,7 @@ def resample_spm( # run the Matlab SPM resampling import matlab as ml - eng.resample_spm_m( + eng.amypad_resample( imrefu, imflou, ml.double(M.tolist()), mask, mean, intrp, which, prefix ) @@ -297,16 +327,18 @@ def resample_spm( if del_out_uncmpr: os.remove(fim) - # > the compressed output naming + # the compressed output naming if fimout: fout = os.path.join(opth, fimout) elif pickname == "ref": fout = os.path.join( - opth, "affine_ref-" + nii.file_parts(imrefu)[1] + fcomment + ".nii.gz", + opth, + "affine_ref-" + nii.file_parts(imrefu)[1] + fcomment + ".nii.gz", ) elif pickname == "flo": fout = os.path.join( - opth, "affine_flo-" + nii.file_parts(imflo)[1] + fcomment + ".nii.gz", + opth, + "affine_flo-" + nii.file_parts(imflo)[1] + fcomment + ".nii.gz", ) # change the file name os.rename(fim + ".gz", fout) @@ -320,3 +352,91 @@ def resample_spm( ) return fout + + +def seg_spm( + f_mri, + spm_path=None, + matlab_eng_name="", + outpath=None, + store_nat_gm=False, + store_nat_wm=False, + store_nat_csf=False, + store_fwd=False, + store_inv=False, + visual=False, +): + """ + Normalisation/segmentation using SPM12. + Args: + f_mri: file path to the T1w MRI file + spm_path(str): SPM path + matlab_eng_name: name of the Python engine for Matlab. + outpath: output folder path for the normalisation file output + store_nat_*: stores native space segmentation output for either + grey matter, white matter or CSF + sotre_fwd/inv: stores forward/inverse normalisation definitions + visual: shows the Matlab window progress + """ + out = {} # output dictionary + # get Matlab engine or use the provided one + eng = ensure_spm(matlab_eng_name) + if not spm_path: + spm_path = spm_dir() + # run SPM normalisation/segmentation + param, invdef, fordef = eng.amypad_seg( + f_mri, + str(spm_path), + float(store_nat_gm), + float(store_nat_wm), + float(store_nat_csf), + float(store_fwd), + float(store_inv), + float(visual), + nargout=3, + ) + if outpath is not None: + create_dir(outpath) + out["param"] = move_files(param, outpath) + out["invdef"] = move_files(invdef, outpath) + out["fordef"] = move_files(fordef, outpath) + # go through tissue types and move them to the output folder + for c in glob_match(r"c\d*", os.path.dirname(param)): + nm = os.path.basename(c)[:2] + out[nm] = move_files(c, outpath) + else: + out["param"] = param + out["invdef"] = invdef + out["fordef"] = fordef + + for c in glob_match(r"c\d*", os.path.dirname(param)): + nm = os.path.basename(c)[:2] + out[nm] = c + return out + + +def normw_spm(f_def, files4norm, matlab_eng_name="", outpath=None): + """ + Write normalisation output to NIfTI files using SPM12. + Args: + f_def: NIfTI file of definitions for non-rigid normalisation + files4norm: list of input NIfTI file paths in format ['file, 1'] + matlab_eng_name: name of the Python engine for Matlab. + outpath: output folder path for the normalisation files + """ + eng = ensure_spm(matlab_eng_name) # get_matlab + eng.amypad_normw(f_def, files4norm) + out = [] # output list + if outpath is not None: + create_dir(outpath) + for f in files4norm: + fpth = f.split(",")[0] + out.append( + move_files( + os.path.join(os.path.dirname(fpth), "w" + os.path.basename(fpth)), + outpath, + ) + ) + else: + out.append("w" + os.path.basename(f.split(",")[0])) + return out diff --git a/spm12/utils.py b/spm12/utils.py index da67894..afbadbb 100644 --- a/spm12/utils.py +++ b/spm12/utils.py @@ -13,7 +13,7 @@ except ImportError: # fix py2.7 from backports.functools_lru_cache import lru_cache -__all__ = ["get_matlab", "ensure_spm"] +__all__ = ["ensure_spm", "get_matlab", "spm_dir"] PATH_M = resource_filename(__name__, "") log = logging.getLogger(__name__) @@ -26,14 +26,19 @@ def get_matlab(name=None): return eng +def spm_dir(cache="~/.spm12", version=12): + cache = path.expanduser(cache) + if str(version) != "12": + raise NotImplementedError + return path.join(cache, "spm12") + + @lru_cache() @wraps(get_matlab) def ensure_spm(name=None, cache="~/.spm12", version=12): eng = get_matlab(name) cache = path.expanduser(cache) - if str(version) != "12": - raise NotImplementedError - addpath = path.join(cache, "spm12") + addpath = spm_dir(cache=cache, version=version) if path.exists(addpath): eng.addpath(addpath) if not eng.exist("spm_jobman"): diff --git a/tests/test_regseg.py b/tests/test_regseg.py index a3865e8..002cc87 100644 --- a/tests/test_regseg.py +++ b/tests/test_regseg.py @@ -34,7 +34,14 @@ def assert_equal_arrays(x, y, nmse_tol=0, denan=True): y: {:.3g}/{:.3g}/{:.3g}({:.3g}) """ ).format( - x.min(), x.mean(), x.max(), x.std(), y.min(), y.mean(), y.max(), y.std(), + x.min(), + x.mean(), + x.max(), + x.std(), + y.min(), + y.mean(), + y.max(), + y.std(), ) )