本文整理汇总了Python中nipy.modalities.fmri.glm.FMRILinearModel.contrast方法的典型用法代码示例。如果您正苦于以下问题:Python FMRILinearModel.contrast方法的具体用法?Python FMRILinearModel.contrast怎么用?Python FMRILinearModel.contrast使用的例子?那么恭喜您, 这里精选的方法代码示例或许可以为您提供帮助。您也可以进一步了解该方法所在类nipy.modalities.fmri.glm.FMRILinearModel
的用法示例。
在下文中一共展示了FMRILinearModel.contrast方法的15个代码示例,这些例子默认根据受欢迎程度排序。您可以为喜欢或者感觉有用的代码点赞,您的评价将有助于系统推荐出更棒的Python代码示例。
示例1: tortoise
# 需要导入模块: from nipy.modalities.fmri.glm import FMRILinearModel [as 别名]
# 或者: from nipy.modalities.fmri.glm.FMRILinearModel import contrast [as 别名]
def tortoise(*args):
print args
print (
'Fitting a "Fixed Effect" GLM for merging LR and RL '
'phase-encoding directions for subject %s ...' % (
subject_data.subject_id))
fmri_glm = FMRILinearModel(subject_data.func,
[design_matrix.matrix
for design_matrix in design_matrices],
mask='compute'
)
fmri_glm.fit(do_scaling=True, model='ar1')
print "... done.\r\n"
# save computed mask
mask_path = os.path.join(subject_data.output_dir, "mask.nii")
print "Saving mask image to %s ..." % mask_path
nibabel.save(fmri_glm.mask, mask_path)
print "... done.\r\n"
z_maps = {}
effects_maps = {}
map_dirs = {}
try:
for contrast_id, contrast_val in contrasts.iteritems():
print "\tcontrast id: %s" % contrast_id
z_map, eff_map = fmri_glm.contrast(
contrast_val,
con_id=contrast_id,
output_z=True,
output_effects=True
)
# store stat maps to disk
for map_type, out_map in zip(['z', 'effects'],
[z_map, eff_map]):
map_dir = os.path.join(
subject_data.output_dir, '%s_maps' % map_type)
map_dirs[map_type] = map_dir
if not os.path.exists(map_dir):
os.makedirs(map_dir)
map_path = os.path.join(
map_dir, '%s_%s.nii' % (map_type, contrast_id))
print "\t\tWriting %s ..." % map_path
nibabel.save(out_map, map_path)
# collect zmaps for contrasts we're interested in
if map_type == 'z':
z_maps[contrast_id] = map_path
if map_type == 'effects':
effects_maps[contrast_id] = map_path
return effects_maps, z_maps, mask_path, map_dirs
except:
return None
示例2: _first_level_glm
# 需要导入模块: from nipy.modalities.fmri.glm import FMRILinearModel [as 别名]
# 或者: from nipy.modalities.fmri.glm.FMRILinearModel import contrast [as 别名]
def _first_level_glm(study_dir, subject_id, model_id,
hrf_model='canonical', drift_model='cosine',
glm_model='ar1', mask='compute', verbose=1):
study_id = os.path.split(study_dir)[1]
subject_dir = os.path.join(study_dir, subject_id)
if verbose > 0:
print '%[email protected]%s: first level glm' % (subject_id, study_id)
tr = get_study_tr(study_dir)
images, n_scans = get_subject_bold_images(subject_dir)
motion = get_subject_motion_per_session(subject_dir)
contrasts = get_task_contrasts(study_dir, subject_dir, model_id)
events = get_subject_events(study_dir, subject_dir)
design_matrices = make_design_matrices(events, n_scans, tr,
hrf_model, drift_model, motion)
glm = FMRILinearModel(images, design_matrices, mask=mask)
glm.fit(do_scaling=True, model=glm_model)
for contrast_id in contrasts:
con_val = []
for session_con, session_dm in zip(contrasts[contrast_id],
design_matrices):
con = np.zeros(session_dm.shape[1])
con[:len(session_con)] = session_con
con_val.append(con)
z_map, t_map, c_map, var_map = glm.contrast(
con_val,
con_id=contrast_id,
output_z=True,
output_stat=True,
output_effects=True,
output_variance=True,)
model_dir = os.path.join(subject_dir, 'model', model_id)
for dtype, img in zip(['z', 't', 'c', 'var'],
[z_map, t_map, c_map, var_map]):
map_dir = os.path.join(model_dir, '%s_maps' % dtype)
if not os.path.exists(map_dir):
os.makedirs(map_dir)
path = os.path.join(
map_dir, '%s.nii.gz' % normalize_name(contrast_id))
nb.save(img, path)
nb.save(glm.mask, os.path.join(model_dir, 'mask.nii.gz'))
示例3: _apply_glm
# 需要导入模块: from nipy.modalities.fmri.glm import FMRILinearModel [as 别名]
# 或者: from nipy.modalities.fmri.glm.FMRILinearModel import contrast [as 别名]
def _apply_glm(out_dir, data, design_matrices,
contrasts, mask='compute', model_id=None, resample=True):
# print out_dir
bold_dir = os.path.join(out_dir, 'fmri')
if not os.path.exists(bold_dir):
os.makedirs(bold_dir)
for i, img in enumerate(data):
if type(img) is str:
img = nb.load(img)
nb.save(img, os.path.join(bold_dir, 'bold_session_%i.nii.gz' % i))
# fit glm
glm = FMRILinearModel(data, design_matrices, mask=mask)
glm.fit(do_scaling=True, model='ar1')
nb.save(glm.mask, os.path.join(out_dir, 'mask.nii.gz'))
if resample:
resample_niimg(os.path.join(out_dir, 'mask.nii.gz'))
stat_maps = {}
for contrast_id in contrasts:
stat_maps[contrast_id] = {}
z_map, t_map, c_map, var_map = glm.contrast(
contrasts[contrast_id],
con_id=contrast_id,
output_z=True,
output_stat=True,
output_effects=True,
output_variance=True,)
for dtype, out_map in zip(['z', 't', 'c', 'variance'],
[z_map, t_map, c_map, var_map]):
map_dir = os.path.join(out_dir, '%s_maps' % dtype)
if not os.path.exists(map_dir):
os.makedirs(map_dir)
if model_id:
map_path = os.path.join(
map_dir, '%s_%s.nii.gz' % (model_id, contrast_id))
else:
map_path = os.path.join(
map_dir, '%s.nii.gz' % contrast_id)
nb.save(out_map, map_path)
if resample:
resample_niimg(map_path)
stat_maps[contrast_id][dtype] = map_path
return stat_maps
示例4: _first_level
# 需要导入模块: from nipy.modalities.fmri.glm import FMRILinearModel [as 别名]
# 或者: from nipy.modalities.fmri.glm.FMRILinearModel import contrast [as 别名]
def _first_level(out_dir, data, design_matrices, contrasts,
glm_model='ar1', mask='compute', verbose=1):
if verbose:
print '%s:' % out_dir
data = check_niimgs(data)
design_matrices = check_design_matrices(design_matrices)
contrasts = check_contrasts(contrasts)
glm = FMRILinearModel(data, design_matrices, mask=mask)
glm.fit(do_scaling=True, model=glm_model)
for i, contrast_id in enumerate(contrasts):
if verbose:
print ' %s/%s - %s ' % (i, len(contrasts), contrast_id)
con_val = []
for session_con, session_dm in zip(contrasts[contrast_id],
design_matrices):
con = np.zeros(session_dm.shape[1])
con[:len(session_con)] = session_con
con_val.append(con)
z_map, t_map, c_map, var_map = glm.contrast(
con_val,
con_id=contrast_id,
output_z=True,
output_stat=True,
output_effects=True,
output_variance=True,)
for dtype, img in zip(['z', 't', 'c', 'var'],
[z_map, t_map, c_map, var_map]):
map_dir = os.path.join(out_dir, '%s_maps' % dtype)
if not os.path.exists(map_dir):
os.makedirs(map_dir)
path = os.path.join(
map_dir, '%s.nii.gz' % remove_special(contrast_id))
nb.save(img, path)
nb.save(glm.mask, os.path.join(out_dir, 'mask.nii.gz'))
示例5: FMRILinearModel
# 需要导入模块: from nipy.modalities.fmri.glm import FMRILinearModel [as 别名]
# 或者: from nipy.modalities.fmri.glm.FMRILinearModel import contrast [as 别名]
for contrast_id in contrasts:
print "\tcontrast id: %s" % contrast_id
# effects maps will be the input to the second level GLM
first_level_image = nibabel.concat_images(
[x[contrast_id] for x in sujects_effects_maps])
# fit 2nd level GLM for given contrast
group_model = FMRILinearModel(first_level_image,
design_matrix, group_mask)
group_model.fit(do_scaling=False, model='ols')
# specify and estimate the contrast
contrast_val = np.array(([[1.]])) # the only possible contrast !
z_map, = group_model.contrast(contrast_val,
con_id='one_sample %s' % contrast_id,
output_z=True)
# save map
map_dir = os.path.join(output_dir, 'z_maps')
if not os.path.exists(map_dir):
os.makedirs(map_dir)
map_path = os.path.join(map_dir, '2nd_level_%s.nii.gz' % (
contrast_id))
print "\t\tWriting %s ..." % map_path
nibabel.save(z_map, map_path)
group_level_z_maps[contrast_id] = map_path
# do stats report
stats_report_filename = os.path.join(output_dir,
示例6: group_one_sample_t_test
# 需要导入模块: from nipy.modalities.fmri.glm import FMRILinearModel [as 别名]
# 或者: from nipy.modalities.fmri.glm.FMRILinearModel import contrast [as 别名]
def group_one_sample_t_test(masks, effects_maps, contrasts, output_dir,
start_time=base_reporter.pretty_time(),
**kwargs):
"""
Runs a one-sample t-test procedure for group analysis. Here, we are
for each experimental condition, only interested refuting the null
hypothesis H0: "The average effect accross the subjects is zero!"
Parameters
----------
masks: list of strings or nibabel image objects
subject masks, one per subject
effects_maps: list of dicts of lists
effects maps from subject-level GLM; each entry is a dictionary;
each entry (indexed by condition id) of this dictionary is the
filename (or correspinding nibabel image object) for the effects
maps for that condition (aka contrast),for that subject
contrasts: dictionary of array_likes
contrasts vectors, indexed by condition id
kwargs: dict_like
parameters can be regular `nipy.labs.viz.plot_map` parameters
(e.g slicer="y") or any parameter we want be reported (e.g
fwhm=[5, 5, 5])
"""
# make output directory
if not os.path.exists(output_dir):
os.makedirs(output_dir)
assert len(masks) == len(effects_maps), (len(masks), len(effects_maps))
# compute group mask
group_mask = nibabel.Nifti1Image(
intersect_masks(masks).astype(np.int8),
(nibabel.load(masks[0]) if isinstance(
masks[0], basestring) else masks[0]).get_affine())
# construct design matrix (only one covariate, namely the "mean effect")
design_matrix = np.ones(len(effects_maps)
)[:, np.newaxis] # only the intercept
group_level_z_maps = {}
group_level_t_maps = {}
for contrast_id in contrasts:
print "\tcontrast id: %s" % contrast_id
# effects maps will be the input to the second level GLM
first_level_image = nibabel.concat_images(
[x[contrast_id] for x in effects_maps])
# fit 2nd level GLM for given contrast
group_model = FMRILinearModel(first_level_image,
design_matrix, group_mask)
group_model.fit(do_scaling=False, model='ols')
# specify and estimate the contrast
contrast_val = np.array(([[1.]])
) # the only possible contrast !
z_map, t_map = group_model.contrast(contrast_val,
con_id='one_sample %s' % contrast_id,
output_z=True,
output_stat=True)
# save map
for map_type, map_img in zip(["z", "t"], [z_map, t_map]):
map_dir = os.path.join(output_dir, '%s_maps' % map_type)
if not os.path.exists(map_dir):
os.makedirs(map_dir)
map_path = os.path.join(map_dir, 'group_level_%s.nii.gz' % (
contrast_id))
print "\t\tWriting %s ..." % map_path
nibabel.save(map_img, map_path)
if map_type == "z":
group_level_z_maps[contrast_id] = map_path
elif map_type == "t":
group_level_z_maps[contrast_id] = map_path
# do stats report
stats_report_filename = os.path.join(
output_dir, "report_stats.html")
generate_subject_stats_report(stats_report_filename, contrasts,
group_level_z_maps, group_mask,
start_time=start_time,
**kwargs)
print "\r\nStatistic report written to %s\r\n" % (
stats_report_filename)
return group_level_z_maps
示例7: do_glm_for_subject
# 需要导入模块: from nipy.modalities.fmri.glm import FMRILinearModel [as 别名]
# 或者: from nipy.modalities.fmri.glm.FMRILinearModel import contrast [as 别名]
def do_glm_for_subject(subject_id, bold_base_folder, trial_base_folder,
output_base_folder):
subject_dir = path(bold_base_folder) / ("sub%03d" % subject_id)
output_dir = (path(output_base_folder) / ("sub%03d" % subject_id) /
"model001")
print output_dir
if not output_dir.exists():
output_dir.makedirs()
anat_file = subject_dir / "highres001.nii"
anat = nb.load(anat_file)
run_ids = range(1, 10)
task_bold_files = [subject_dir.glob("task001_run%03d/rbold*.nii"
% rid)[0]
for rid in run_ids]
task_mvt_files = [subject_dir.glob("task001_run%03d/rp_bold*.txt" %
rid)[0]
for rid in run_ids]
trial_files = [(path(trial_base_folder) / ("Sub%02d" % subject_id) /
"BOLD" / "Trials" / ("run_%02d_spmdef.txt" % rid))
for rid in range(1, 10)]
stats_start_time = pretty_time()
paradigms = []
design_matrices = []
n_scans = []
all_frametimes = []
list_of_contrast_dicts = [] # one dict per run
for bold_file, mvt_file, trial_file in zip(task_bold_files,
task_mvt_files,
trial_files):
_n_scans = nb.load(bold_file).shape[-1]
n_scans.append(_n_scans)
paradigm = make_paradigm(trial_file)
paradigms.append(paradigm)
movements = np.loadtxt(mvt_file)
tr = 2.
drift_model = "Cosine"
hrf_model = "Canonical With Derivative"
hfcut = 128.
frametimes = np.linspace(0, (_n_scans - 1) * tr, _n_scans)
design_matrix = make_dmtx(
frametimes,
paradigm,
hrf_model=hrf_model,
drift_model=drift_model,
hfcut=hfcut,
add_regs=movements,
add_reg_names=[
"Tx", "Ty", "Tz", "R1", "R2", "R3"])
design_matrices.append(design_matrix)
all_frametimes.append(frametimes)
# specify contrasts
contrasts = {}
n_columns = len(design_matrix.names)
for i in xrange(paradigm.n_conditions):
contrasts['%s' % design_matrix.names[2 * i]] = np.eye(
n_columns)[2 * i]
# more interesting contrasts"""
contrasts['Famous-Unfamiliar'] = contrasts[
'Famous'] - contrasts['Unfamiliar']
contrasts['Unfamiliar-Famous'] = -contrasts['Famous-Unfamiliar']
contrasts['Famous-Scrambled'] = contrasts[
'Famous'] - contrasts['Scrambled']
contrasts['Scrambled-Famous'] = -contrasts['Famous-Scrambled']
contrasts['Unfamiliar-Scrambled'] = contrasts[
'Unfamiliar'] - contrasts['Scrambled']
contrasts['Scrambled-Unfamiliar'] = -contrasts['Unfamiliar-Scrambled']
list_of_contrast_dicts.append(contrasts)
# importat maps
z_maps = {}
effects_maps = {}
fmri_glm = FMRILinearModel(task_bold_files,
[dm.matrix for dm in design_matrices],
mask="compute")
fmri_glm.fit(do_scaling=True, model="ar1")
# replicate contrasts across runs
contrasts = dict((cid, [contrasts[cid]
for contrasts in list_of_contrast_dicts])
for cid, cval in contrasts.iteritems())
# compute effects
for contrast_id, contrast_val in contrasts.iteritems():
print "\tcontrast id: %s" % contrast_id
#.........这里部分代码省略.........
示例8: do_subject_glm
# 需要导入模块: from nipy.modalities.fmri.glm import FMRILinearModel [as 别名]
# 或者: from nipy.modalities.fmri.glm.FMRILinearModel import contrast [as 别名]
def do_subject_glm(subject_data):
"""FE analysis for a single subject."""
subject_id = subject_data['subject_id']
output_dir = subject_data["output_dir"]
func_files = subject_data['func']
anat = subject_data['anat']
onset_files = subject_data['onset']
# subject_id = os.path.basename(subject_dir)
# subject_output_dir = os.path.join(output_dir, subject_id)
mem = Memory(os.path.join(output_dir, "cache"))
if not os.path.exists(output_dir):
os.makedirs(output_dir)
# glob files: anat, session func files, session onset files
# anat = glob.glob(os.path.join(subject_dir, anat_wildcard))
# assert len(anat) == 1
# anat = anat[0]
# onset_files = sorted([glob.glob(os.path.join(subject_dir, session))[0]
# for session in session_onset_wildcards])
# func_files = sorted([sorted(glob.glob(os.path.join(subject_dir, session)))
# for session in session_func_wildcards])
### Preprocess data #######################################################
if 0:
subject_data = mem.cache(do_subject_preproc)(
dict(func=func_files, anat=anat, output_dir=output_dir))
func_files = subject_data['func']
anat = subject_data['anat']
# reslice func images
func_files = [mem.cache(reslice_vols)(
sess_func,
target_affine=nibabel.load(sess_func[0]).get_affine())
for sess_func in func_files]
### GLM: loop on (session_bold, onse_file) pairs over the various sessions
design_matrices = []
for session, (func_file, onset_file) in enumerate(zip(func_files,
onset_files)):
if isinstance(func_file, str):
bold = nibabel.load(func_file)
else:
if len(func_file) == 1:
func_file = func_file[0]
bold = nibabel.load(func_file)
assert len(bold.shape) == 4
n_scans = bold.shape[-1]
del bold
else:
n_scans = len(func_file)
frametimes = np.linspace(0, (n_scans - 1) * tr, n_scans)
conditions, onsets, durations, amplitudes = parse_onset_file(
onset_file)
onsets *= tr
durations *= tr
paradigm = BlockParadigm(con_id=conditions, onset=onsets,
duration=durations, amplitude=amplitudes)
design_matrices.append(make_dmtx(frametimes,
paradigm, hrf_model=hrf_model,
drift_model=drift_model,
hfcut=hfcut))
# specify contrasts
n_columns = len(design_matrices[0].names)
contrasts = {}
for i in xrange(paradigm.n_conditions):
contrasts['%s' % design_matrices[0].names[2 * i]
] = np.eye(n_columns)[2 * i]
# more interesting contrasts
contrasts['faces-scrambled'] = contrasts['faces'
] - contrasts['scrambled']
contrasts['scrambled-faces'] = -contrasts['faces-scrambled']
contrasts['effects_of_interest'] = contrasts['faces'
] + contrasts['scrambled']
# effects of interest F-test
diff_contrasts = []
for i in xrange(paradigm.n_conditions - 1):
a = contrasts[design_matrices[0].names[2 * i]]
b = contrasts[design_matrices[0].names[2 * (i + 1)]]
diff_contrasts.append(a - b)
contrasts["diff"] = diff_contrasts
# fit GLM
print 'Fitting a GLM (this takes time)...'
fmri_glm = FMRILinearModel([nibabel.concat_images(sess_func,
check_affines=False)
for sess_func in func_files],
[design_matrix.matrix
for design_matrix in design_matrices],
mask='compute'
)
fmri_glm.fit(do_scaling=True, model='ar1')
# save computed mask
mask_path = os.path.join(output_dir, "mask.nii.gz")
print "Saving mask image %s" % mask_path
nibabel.save(fmri_glm.mask, mask_path)
#.........这里部分代码省略.........
示例9: GLM
# 需要导入模块: from nipy.modalities.fmri.glm import FMRILinearModel [as 别名]
# 或者: from nipy.modalities.fmri.glm.FMRILinearModel import contrast [as 别名]
print 'Fitting a GLM (this takes time)...'
fmri_glm = FMRILinearModel(data_path, design_matrix.matrix,
mask='compute')
fmri_glm.fit(do_scaling=True, model='ar1')
#########################################
# Estimate the contrasts
#########################################
print 'Computing contrasts...'
for index, (contrast_id, contrast_val) in enumerate(contrasts.iteritems()):
print ' Contrast % 2i out of %i: %s' % (
index + 1, len(contrasts), contrast_id)
# save the z_image
image_path = path.join(write_dir, '%s_z_map.nii' % contrast_id)
z_map, = fmri_glm.contrast(contrast_val, con_id=contrast_id, output_z=True)
save(z_map, image_path)
# Create snapshots of the contrasts
vmax = max(- z_map.get_data().min(), z_map.get_data().max())
if index > 0:
plt.clf()
plot_map(z_map.get_data(), z_map.get_affine(),
cmap=cm.cold_hot,
vmin=- vmax,
vmax=vmax,
anat=None,
figure=10,
threshold=2.5)
plt.savefig(path.join(write_dir, '%s_z_map.png' % contrast_id))
示例10: FMRILinearModel
# 需要导入模块: from nipy.modalities.fmri.glm import FMRILinearModel [as 别名]
# 或者: from nipy.modalities.fmri.glm.FMRILinearModel import contrast [as 别名]
# Input files
fmri_files = [example_data.get_filename('fiac', 'fiac0', run)
for run in ['run1.nii.gz', 'run2.nii.gz']]
design_files = [example_data.get_filename('fiac', 'fiac0', run)
for run in ['run1_design.npz', 'run2_design.npz']]
mask_file = example_data.get_filename('fiac', 'fiac0', 'mask.nii.gz')
# Load all the data
multi_session_model = FMRILinearModel(fmri_files, design_files, mask_file)
# GLM fitting
multi_session_model.fit(do_scaling=True, model='ar1')
# Compute the required contrast
print('Computing test contrast image...')
n_regressors = [np.load(f)['X'].shape[1] for f in design_files]
con = [np.hstack((cvect, np.zeros(nr - len(cvect)))) for nr in n_regressors]
z_map, = multi_session_model.contrast(con)
# Show Z-map image
mean_map = multi_session_model.means[0]
plot_map(z_map.get_data(),
z_map.get_affine(),
anat=mean_map.get_data(),
anat_affine=mean_map.get_affine(),
cmap=cm.cold_hot,
threshold=2.5,
black_bg=True)
plt.show()
示例11: Nifti1Image
# 需要导入模块: from nipy.modalities.fmri.glm import FMRILinearModel [as 别名]
# 或者: from nipy.modalities.fmri.glm.FMRILinearModel import contrast [as 别名]
grp_mask = Nifti1Image(intersect_masks(mask_images).astype(np.int8),
load(mask_images[0]).get_affine())
# concatenate the individual images
first_level_image = concat_images(betas)
# set the model
design_matrix = np.ones(len(betas))[:, np.newaxis] # only the intercept
grp_model = FMRILinearModel(first_level_image, design_matrix, grp_mask)
# GLM fitting using ordinary least_squares
grp_model.fit(do_scaling=False, model='ols')
# specify and estimate the contrast
contrast_val = np.array(([[1]])) # the only possible contrast !
z_map, = grp_model.contrast(contrast_val, con_id='one_sample', output_z=True)
# write the results
save(z_map, path.join(write_dir, 'one_sample_z_map.nii'))
# look at the result
vmax = max(- z_map.get_data().min(), z_map.get_data().max())
vmin = - vmax
plot_map(z_map.get_data(), z_map.get_affine(),
cmap=cm.cold_hot,
vmin=vmin,
vmax=vmax,
threshold=3.,
black_bg=True)
plt.savefig(path.join(write_dir, '%s_z_map.png' % 'one_sample'))
plt.show()
示例12: first_level
# 需要导入模块: from nipy.modalities.fmri.glm import FMRILinearModel [as 别名]
# 或者: from nipy.modalities.fmri.glm.FMRILinearModel import contrast [as 别名]
def first_level(subject_dic):
# experimental paradigm meta-params
stats_start_time = time.ctime()
tr = 2.4
drift_model = 'blank'
hrf_model = 'canonical' # hemodynamic reponse function
hfcut = 128.
n_scans = 128
# make design matrices
mask_images = []
design_matrices = []
fmri_files = subject_dic['func']
for x in xrange(len(fmri_files)):
paradigm = paradigm_contrasts.localizer_paradigm()
# build design matrix
frametimes = np.linspace(0, (n_scans - 1) * tr, n_scans)
design_matrix = make_dmtx(
frametimes,
paradigm, hrf_model=hrf_model,
drift_model=drift_model, hfcut=hfcut,
)
design_matrices.append(design_matrix)
# Specify contrasts
contrasts = paradigm_contrasts.localizer_contrasts(design_matrix)
#create output directory
subject_session_output_dir = os.path.join(subject_dic['output_dir'],
'res_stats')
if not os.path.exists(subject_session_output_dir):
os.makedirs(subject_session_output_dir)
# Fit GLM
print 'Fitting a GLM (this takes time)...'
fmri_glm = FMRILinearModel(fmri_files,
[design_matrix.matrix
for design_matrix in design_matrices],
mask='compute'
)
fmri_glm.fit(do_scaling=True, model='ar1')
# save computed mask
mask_path = os.path.join(subject_session_output_dir,
"mask.nii.gz")
print "Saving mask image %s" % mask_path
nibabel.save(fmri_glm.mask, mask_path)
mask_images.append(mask_path)
# compute contrasts
z_maps = {}
effects_maps = {}
for contrast_id, contrast_val in contrasts.iteritems():
print "\tcontrast id: %s" % contrast_id
z_map, t_map, effects_map, var_map = fmri_glm.contrast(
[contrast_val] * 1,
con_id=contrast_id,
output_z=True,
output_stat=True,
output_effects=True,
output_variance=True
)
# store stat maps to disk
for map_type, out_map in zip(['z', 't', 'effects', 'variance'],
[z_map, t_map, effects_map, var_map]):
map_dir = os.path.join(
subject_session_output_dir, '%s_maps' % map_type)
if not os.path.exists(map_dir):
os.makedirs(map_dir)
map_path = os.path.join(
map_dir, '%s%s.nii.gz' %(subject_dic['subject_id'], contrast_id))
print "\t\tWriting %s ..." % map_path
nibabel.save(out_map, map_path)
# collect zmaps for contrasts we're interested in
if map_type == 'z':
z_maps[contrast_id] = map_path
if map_type == 'effects':
effects_maps[contrast_id] = map_path
# do stats report
anat_img = nibabel.load(subject_dic['anat'])
stats_report_filename = os.path.join(subject_session_output_dir,
"report_stats.html")
generate_subject_stats_report(
stats_report_filename,
contrasts,
z_maps,
fmri_glm.mask,
threshold=2.3,
cluster_th=15,
anat=anat_img,
anat_affine=anat_img.get_affine(),
design_matrices=design_matrix,
subject_id="sub001",
#.........这里部分代码省略.........
示例13: zip
# 需要导入模块: from nipy.modalities.fmri.glm import FMRILinearModel [as 别名]
# 或者: from nipy.modalities.fmri.glm.FMRILinearModel import contrast [as 别名]
mask='compute')
fmri_glm.fit(do_scaling=True, model='ar1')
# save computed mask
mask_path = os.path.join(subject_data.output_dir, "mask.nii.gz")
print "Saving mask image %s" % mask_path
nibabel.save(fmri_glm.mask, mask_path)
mask_images.append(mask_path)
# compute contrast maps
z_maps = {}
effects_maps = {}
for contrast_id, contrast_val in contrasts.iteritems():
print "\tcontrast id: %s" % contrast_id
z_map, t_map, effects_map, var_map = fmri_glm.contrast(
[contrast_val] * 2, con_id=contrast_id, output_z=True,
output_stat=True, output_effects=True, output_variance=True)
for map_type, out_map in zip(['z', 't', 'effects', 'variance'],
[z_map, t_map, effects_map, var_map]):
map_dir = os.path.join(
subject_data.output_dir, '%s_maps' % map_type)
if not os.path.exists(map_dir):
os.makedirs(map_dir)
map_path = os.path.join(
map_dir, '%s.nii.gz' % contrast_id)
print "\t\tWriting %s ..." % map_path
nibabel.save(out_map, map_path)
if map_type == 'z':
z_maps[contrast_id] = map_path
if map_type == 'effects':
effects_maps[contrast_id] = map_path
示例14: do_subject_glm
# 需要导入模块: from nipy.modalities.fmri.glm import FMRILinearModel [as 别名]
# 或者: from nipy.modalities.fmri.glm.FMRILinearModel import contrast [as 别名]
def do_subject_glm(subject_data):
"""FE analysis for a single subject."""
subject_id = subject_data['subject_id']
output_dir = subject_data["output_dir"]
func_files = subject_data['func']
anat = subject_data['anat']
onset_files = subject_data['onset']
tr = subject_data['TR']
time_units = subject_data['time_units'].lower()
assert time_units in ["seconds", "tr", "milliseconds"]
drift_model = subject_data['drift_model']
hrf_model = subject_data["hrf_model"]
hfcut = subject_data["hfcut"]
mem = Memory(os.path.join(output_dir, "cache"))
if not os.path.exists(output_dir):
os.makedirs(output_dir)
if 0:
subject_data = mem.cache(do_subject_preproc)(
dict(func=func_files, anat=anat, output_dir=output_dir))
func_files = subject_data['func']
anat = subject_data['anat']
# reslice func images
func_files = [mem.cache(reslice_vols)(
sess_func, target_affine=nibabel.load(sess_func[0]).get_affine())
for sess_func in func_files]
### GLM: loop on (session_bold, onse_file) pairs over the various sessions
design_matrices = []
for func_file, onset_file in zip(func_files, onset_files):
if isinstance(func_file, str):
bold = nibabel.load(func_file)
else:
if len(func_file) == 1:
func_file = func_file[0]
bold = nibabel.load(func_file)
assert len(bold.shape) == 4
n_scans = bold.shape[-1]
del bold
else:
n_scans = len(func_file)
frametimes = np.linspace(0, (n_scans - 1) * tr, n_scans)
conditions, onsets, durations, amplitudes = parse_onset_file(
onset_file)
if time_units == "tr":
onsets *= tr
durations *= tr
elif time_units in ["milliseconds"]:
onsets *= 1e-3
durations *= 1e-3
paradigm = BlockParadigm(con_id=conditions, onset=onsets,
duration=durations, amplitude=amplitudes)
design_matrices.append(make_dmtx(
frametimes,
paradigm, hrf_model=hrf_model,
drift_model=drift_model, hfcut=hfcut))
# specify contrasts
n_columns = len(design_matrices[0].names)
contrasts = {}
for i in range(paradigm.n_conditions):
contrasts['%s' % design_matrices[0].names[2 * i]
] = np.eye(n_columns)[2 * i]
# effects of interest F-test
diff_contrasts = []
for i in range(paradigm.n_conditions - 1):
a = contrasts[design_matrices[0].names[2 * i]]
b = contrasts[design_matrices[0].names[2 * (i + 1)]]
diff_contrasts.append(a - b)
contrasts["diff"] = diff_contrasts
# fit GLM
print('Fitting a GLM (this takes time)...')
fmri_glm = FMRILinearModel([nibabel.concat_images(sess_func,
check_affines=False)
for sess_func in func_files],
[design_matrix.matrix
for design_matrix in design_matrices],
mask='compute'
)
fmri_glm.fit(do_scaling=True, model='ar1')
# save computed mask
mask_path = os.path.join(output_dir, "mask.nii.gz")
print("Saving mask image %s" % mask_path)
nibabel.save(fmri_glm.mask, mask_path)
# compute contrasts
z_maps = {}
effects_maps = {}
for contrast_id, contrast_val in contrasts.items():
print("\tcontrast id: %s" % contrast_id)
if np.ndim(contrast_val) > 1:
contrast_type = "t"
else:
contrast_type = "F"
z_map, t_map, effects_map, var_map = fmri_glm.contrast(
#.........这里部分代码省略.........
示例15: range
# 需要导入模块: from nipy.modalities.fmri.glm import FMRILinearModel [as 别名]
# 或者: from nipy.modalities.fmri.glm.FMRILinearModel import contrast [as 别名]
for i in range(paradigm.n_conditions):
contrasts['%s' % design_matrix.names[2 * i]] = np.eye(n_columns)[2 * i]
# Our contrast of interest
reading_vs_visual = contrasts["phrasevideo"] - contrasts["damier_H"]
########################################
# Perform a GLM analysis on H1
########################################
fmri_glm = FMRILinearModel(fmri_data,
design_matrix.matrix, mask='compute')
fmri_glm.fit(do_scaling=True, model='ar1')
# Estimate the contrast
z_map, = fmri_glm.contrast(reading_vs_visual, output_z=True)
# Plot the contrast
vmax = max(-z_map.get_data().min(), z_map.get_data().max())
plot_map(z_map.get_data(), z_map.get_affine(),
cmap=cm.cold_hot, vmin=-vmax, vmax=vmax,
slicer='z', black_bg=True, threshold=2.5,
title='Reading vs visual')
# Count all the clusters for |Z| > 2.5
Z = z_map.get_data()
from scipy import ndimage
cluster_map, n_clusters = ndimage.label(np.abs(Z) > 2.5)
cluster_sizes = np.bincount(cluster_map.ravel())[1:]
print "Cluster sizes:"