Skip to content

Commit

Permalink
Merge pull request #55 from pnlbwh/fix-force
Browse files Browse the repository at this point in the history 
Fix force
  • Loading branch information
@tashrifbillah
tashrifbillah authored Feb 24, 2020
2 parents bb8e3c7 + 8106997 commit 3187882
Show file tree
Hide file tree
Showing 12 changed files with 265 additions and 252 deletions.
45 changes: 15 additions & 30 deletions README.md
View file
Original file line number Diff line number Diff line change
Expand Up @@ -347,42 +347,27 @@ NOTE: running the above test should take roughly an hour.
`./pipeline_test.sh` will download test data*, and run the whole processing pipeline on them.
If the test is successful and complete, you should see the following output on the command line.

CONNECTOM mean FA: 0.6787134267438087
PRISMA mean FA before harmonization: 0.8016994786111754
PRISMA mean FA after harmonization: 0.6335818260201815
02/24/20 11:37,mean meanFA,std meanFA
CONNECTOM,0.6799002461759734,0.007886013534889139
PRISMA_before,0.8016994786111754,0.010304926419301603
PRISMA_after,0.6763777736686049,0.010977445592601354


In detail:

CONNECTOM site:
mean FA over IIT_mean_FA_skeleton.nii.gz for all cases:
dwi_A_connectom_st_b1200_resampled.nii.gz 0.6664273368527192
dwi_B_connectom_st_b1200_resampled.nii.gz 0.6880572316659973
dwi_C_connectom_st_b1200_resampled.nii.gz 0.6816557117127096

mean meanFA: 0.6787134267438087
std meanFA: 0.00907215035284852

PRISMA site before harmonization:
mean FA over IIT_mean_FA_skeleton.nii.gz for all cases:
dwi_A_prisma_st_b1200.nii.gz 0.7879066657395326
dwi_B_prisma_st_b1200.nii.gz 0.8126709307616746
dwi_C_prisma_st_b1200.nii.gz 0.8045208393323189

mean meanFA: 0.8016994786111754
std meanFA: 0.010304926419301603
The above output are saved in `template/meanFAstat.csv` file. Detailed statistics are saved in:

connectom_prisma/template/CONNECTOM_stat.csv
connectom_prisma/template/PRISMA_before_stat.csv
connectom_prisma/template/PRISMA_after_stat.csv

PRISMA site after harmonization:
mean FA over IIT_mean_FA_skeleton.nii.gz for all cases:
harmonized_dwi_A_prisma_st_b1200_resampled.nii.gz 0.6240905840733715
harmonized_dwi_B_prisma_st_b1200_resampled.nii.gz 0.6487062683838898
harmonized_dwi_C_prisma_st_b1200_resampled.nii.gz 0.6279486256032834

mean meanFA: 0.6335818260201815
std meanFA: 0.010809954940438949
In addition, you should check out the plot `template/meanFAstat_ebarplot.png` that demonstrates quality of dMRIharmonization:

![](doc/meanFAstat_ebarplot.png)


As you see in the above, after harmonization, target site (PRISMA) meanFA came closer to that of reference site (CONNECTOM) FA.

\* If there is any problem downloading test data, try manually downloading and unzipping it to `lib/tests/` folder.
**NOTE** If there is any problem downloading test data, try manually downloading and unzipping it to `lib/tests/` folder.


## 2. unittest
Expand Down
Binary file added doc/meanFAstat_ebarplot.png
View file
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
88 changes: 44 additions & 44 deletions lib/buildTemplate.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -20,9 +20,9 @@
from util import *

eps= 2.2204e-16
SCRIPTDIR= os.path.dirname(__file__)
config = configparser.ConfigParser()
config.read(f'/tmp/harm_config_{os.getpid()}.ini')
SCRIPTDIR= dirname(__file__)
config = ConfigParser()
config.read(f'/tmp/harm_config_{getpid()}.ini')
N_shm = int(config['DEFAULT']['N_shm'])
N_proc = int(config['DEFAULT']['N_proc'])
diffusionMeasures= [x for x in config['DEFAULT']['diffusionMeasures'].split(',')]
Expand All @@ -41,50 +41,50 @@ def applyXform(inImg, refImg, warp, trans, outImg):

def warp_bands(imgPath, maskPath, templatePath):

prefix= os.path.basename(imgPath).split('.nii')[0]
directory= os.path.dirname(imgPath)
warp = glob(os.path.join(templatePath, prefix + f'_FA*[!Inverse]Warp.nii.gz'))
trans = glob(os.path.join(templatePath, prefix + f'_FA*GenericAffine.mat'))
prefix= basename(imgPath).split('.nii')[0]
directory= dirname(imgPath)
warp = glob(pjoin(templatePath, prefix + f'_FA*[!Inverse]Warp.nii.gz'))
trans = glob(pjoin(templatePath, prefix + f'_FA*GenericAffine.mat'))

# warping the mask
applyXform(maskPath,
os.path.join(templatePath, 'template0.nii.gz'),
pjoin(templatePath, 'template0.nii.gz'),
warp, trans,
os.path.join(templatePath, os.path.abspath(maskPath).split('.nii')[0]+ 'Warped.nii.gz'))
pjoin(templatePath, abspath(maskPath).split('.nii')[0]+ 'Warped.nii.gz'))


# warping the rish features
for i in range(0, N_shm+1, 2):
applyXform(os.path.join(directory, 'harm', f'{prefix}_L{i}.nii.gz'),
os.path.join(templatePath, 'template0.nii.gz'),
applyXform(pjoin(directory, 'harm', f'{prefix}_L{i}.nii.gz'),
pjoin(templatePath, 'template0.nii.gz'),
warp, trans,
os.path.join(templatePath, f'{prefix}_WarpedL{i}.nii.gz'))
pjoin(templatePath, f'{prefix}_WarpedL{i}.nii.gz'))


# warping the diffusion measures
for dm in diffusionMeasures:
applyXform(os.path.join(directory, 'dti', f'{prefix}_{dm}.nii.gz'),
os.path.join(templatePath, 'template0.nii.gz'),
applyXform(pjoin(directory, 'dti', f'{prefix}_{dm}.nii.gz'),
pjoin(templatePath, 'template0.nii.gz'),
warp, trans,
os.path.join(templatePath, f'{prefix}_Warped{dm}.nii.gz'))
pjoin(templatePath, f'{prefix}_Warped{dm}.nii.gz'))


def createAntsCaselist(imgs, file):

with open(file,'w') as f:
for imgPath in imgs:
prefix= os.path.basename(imgPath).split('.nii')[0]
directory= os.path.dirname(imgPath)
prefix= basename(imgPath).split('.nii')[0]
directory= dirname(imgPath)

FA= os.path.join(directory,'dti', f'{prefix}_FA.nii.gz')
L0= os.path.join(directory,'harm', f'{prefix}_L0.nii.gz')
FA= pjoin(directory,'dti', f'{prefix}_FA.nii.gz')
L0= pjoin(directory,'harm', f'{prefix}_L0.nii.gz')
f.write(f'{FA},{L0}\n')


def antsMult(caselist, outPrefix):

N_core=os.getenv('TEMPLATE_CONSTRUCT_CORES')
check_call((' ').join([os.path.join(SCRIPTDIR, 'antsMultivariateTemplateConstruction2_fixed_random_seed.sh'),
N_core=getenv('TEMPLATE_CONSTRUCT_CORES')
check_call((' ').join([pjoin(SCRIPTDIR, 'antsMultivariateTemplateConstruction2_fixed_random_seed.sh'),
'-d', '3',
'-g', '0.2',
'-k', '2',
Expand All @@ -101,25 +101,25 @@ def dti_stat(siteName, imgs, masks, templatePath, templateHdr):

maskData = []
for maskPath in masks:
maskData.append(load_nifti(os.path.join(templatePath, os.path.abspath(maskPath).split('.nii')[0] + 'Warped.nii.gz'))[0])
maskData.append(load_nifti(pjoin(templatePath, abspath(maskPath).split('.nii')[0] + 'Warped.nii.gz'))[0])


morphed_mask= binary_opening(np.mean(maskData, axis= 0)>0.5, structure= generate_binary_structure(3,1))*1
morphed_mask_name= os.path.join(templatePath, f'{siteName}_Mask.nii.gz')
morphed_mask_name= pjoin(templatePath, f'{siteName}_Mask.nii.gz')
templateAffine = templateHdr.get_best_affine()
save_nifti(morphed_mask_name, morphed_mask.astype('uint8'), templateAffine, templateHdr)


for dm in diffusionMeasures:
imgData= []
for imgPath in imgs:
prefix = os.path.basename(imgPath).split('.nii')[0]
imgData.append(load_nifti(os.path.join(templatePath, f'{prefix}_Warped{dm}.nii.gz'))[0])
prefix = basename(imgPath).split('.nii')[0]
imgData.append(load_nifti(pjoin(templatePath, f'{prefix}_Warped{dm}.nii.gz'))[0])

save_nifti(os.path.join(templatePath, f'Mean_{siteName}_{dm}.nii.gz'),
save_nifti(pjoin(templatePath, f'Mean_{siteName}_{dm}.nii.gz'),
np.mean(imgData, axis= 0), templateAffine, templateHdr)

save_nifti(os.path.join(templatePath, f'Std_{siteName}_{dm}.nii.gz'),
save_nifti(pjoin(templatePath, f'Std_{siteName}_{dm}.nii.gz'),
np.std(imgData, axis= 0), templateAffine, templateHdr)

return morphed_mask_name
Expand All @@ -130,14 +130,14 @@ def rish_stat(siteName, imgs, templatePath, templateHdr):
for i in range(0, N_shm+1, 2):
imgData= []
for imgPath in imgs:
prefix = os.path.basename(imgPath).split('.nii')[0]
imgData.append(load_nifti(os.path.join(templatePath, f'{prefix}_WarpedL{i}.nii.gz'))[0])
prefix = basename(imgPath).split('.nii')[0]
imgData.append(load_nifti(pjoin(templatePath, f'{prefix}_WarpedL{i}.nii.gz'))[0])

templateAffine= templateHdr.get_best_affine()
save_nifti(os.path.join(templatePath, f'Mean_{siteName}_L{i}.nii.gz'),
save_nifti(pjoin(templatePath, f'Mean_{siteName}_L{i}.nii.gz'),
np.mean(imgData, axis= 0), templateAffine, templateHdr)

save_nifti(os.path.join(templatePath, f'Std_{siteName}_L{i}.nii.gz'),
save_nifti(pjoin(templatePath, f'Std_{siteName}_L{i}.nii.gz'),
np.std(imgData, axis= 0), templateAffine, templateHdr)


Expand All @@ -148,14 +148,14 @@ def template_masking(refMaskPath, targetMaskPath, templatePath, siteName):

templateMask= applymask(ref.get_data(), target.get_data())

save_nifti(os.path.join(templatePath, 'templateMask.nii.gz'), templateMask.astype('uint8'), ref.affine, ref.header)
save_nifti(pjoin(templatePath, 'templateMask.nii.gz'), templateMask.astype('uint8'), ref.affine, ref.header)

for dm in diffusionMeasures:
fileName= os.path.join(templatePath, f'Mean_{siteName}_{dm}.nii.gz')
fileName= pjoin(templatePath, f'Mean_{siteName}_{dm}.nii.gz')
img= load(fileName)
save_nifti(fileName, applymask(img.get_data(), templateMask), img.affine, img.header)

fileName= os.path.join(templatePath, f'Std_{siteName}_{dm}.nii.gz')
fileName= pjoin(templatePath, f'Std_{siteName}_{dm}.nii.gz')
img= load(fileName)
save_nifti(fileName, applymask(img.get_data(), templateMask), img.affine, img.header)

Expand Down Expand Up @@ -212,11 +212,11 @@ def difference_calc(refSite, targetSite, refImgs, targetImgs,
if travelHeads:
print('Using travelHeads for computing templates of',dm)
for refImg, targetImg in zip(refImgs, targetImgs):
prefix = os.path.basename(refImg).split('.nii')[0]
ref= load_nifti(os.path.join(templatePath, f'{prefix}_Warped{dm}.nii.gz'))[0]
prefix = basename(refImg).split('.nii')[0]
ref= load_nifti(pjoin(templatePath, f'{prefix}_Warped{dm}.nii.gz'))[0]

prefix = os.path.basename(targetImg).split('.nii')[0]
target= load_nifti(os.path.join(templatePath, f'{prefix}_Warped{dm}.nii.gz'))[0]
prefix = basename(targetImg).split('.nii')[0]
target= load_nifti(pjoin(templatePath, f'{prefix}_Warped{dm}.nii.gz'))[0]

temp= stat_calc(ref, target, mask)
delta.append(temp[0])
Expand All @@ -225,10 +225,10 @@ def difference_calc(refSite, targetSite, refImgs, targetImgs,
scale.append(temp[3])

else:
fileName= os.path.join(templatePath, f'Mean_{refSite}_{dm}.nii.gz')
fileName= pjoin(templatePath, f'Mean_{refSite}_{dm}.nii.gz')
ref= load_nifti(fileName)[0]

fileName= os.path.join(templatePath, f'Mean_{targetSite}_{dm}.nii.gz')
fileName= pjoin(templatePath, f'Mean_{targetSite}_{dm}.nii.gz')
target= load_nifti(fileName)[0]

temp= stat_calc(ref, target, mask)
Expand All @@ -237,17 +237,17 @@ def difference_calc(refSite, targetSite, refImgs, targetImgs,
per_diff_smooth.append(temp[2])
scale.append(temp[3])

save_nifti(os.path.join(templatePath, f'Delta_{dm}.nii.gz'),
save_nifti(pjoin(templatePath, f'Delta_{dm}.nii.gz'),
np.mean(delta, axis= 0), templateAffine, templateHdr)

save_nifti(os.path.join(templatePath, f'PercentageDiff_{dm}.nii.gz'),
save_nifti(pjoin(templatePath, f'PercentageDiff_{dm}.nii.gz'),
np.mean(per_diff, axis= 0), templateAffine, templateHdr)

save_nifti(os.path.join(templatePath, f'PercentageDiff_{dm}smooth.nii.gz'),
save_nifti(pjoin(templatePath, f'PercentageDiff_{dm}smooth.nii.gz'),
np.mean(per_diff_smooth, axis= 0), templateAffine, templateHdr)

if 'L' in dm:
save_nifti(os.path.join(templatePath, f'Scale_{dm}.nii.gz'),
save_nifti(pjoin(templatePath, f'Scale_{dm}.nii.gz'),
np.sqrt(np.mean(scale, axis= 0)), templateAffine, templateHdr)


Expand Down
Loading

0 comments on commit 3187882

Please sign in to comment.