AFNI Message Board

Dear AFNI users-

We are very pleased to announce that the new AFNI Message Board framework is up! Please join us at:

https://discuss.afni.nimh.nih.gov

Existing user accounts have been migrated, so returning users can login by requesting a password reset. New users can create accounts, as well, through a standard account creation process. Please note that these setup emails might initially go to spam folders (esp. for NIH users!), so please check those locations in the beginning.

The current Message Board discussion threads have been migrated to the new framework. The current Message Board will remain visible, but read-only, for a little while.

Sincerely, AFNI HQ

History of AFNI updates  

|
July 10, 2012 03:11PM
Andrew,

I wanted to do this a while ago because I noticed that 3dREMLfit seemed to be numerically unstable in the presence of zero columns, and I eventually ended up writing a little python script to strip zero columns out of the X-matrix. Since it was a quick-and-dirty fix, you'll probably need to clean it up a bit before it works (for instance, I hard coded in the path to the AFNI binaries). To use it, you'll want to run 3dDeconvolve -x1D_stop, then run 1dzerostrip.py on the something.xmat.1D file it produces, then run 3dREMLfit with the modified X-matrix. Part of the trick is that 1dtool.py strips out some of the metadata that 3dREMLfit needs to run, and this script preserves that information (but still strips out some sub-TR information about the stimuli, if that's important to you).

If you're using 3dDeconvolve instead of 3dREMLfit, you can do the same thing except feed columns of the modified X-matrix back into 3dDeconvolve as -stim_files, but that seems pretty unwieldy and it might just be easier to complain to the powers-that-be that 3dDeconvolve should be able to read .xmat.1D files. My personal preference with 3dDeconvolve would be to strip out zero columns manually on one subject and compare that output to the output from -allzero_OK. If they're the same, save yourself the bother and just use -allzero_OK.

#!/usr/bin/python

# 1dzerostrip.py
# plagiarized with 'improvements' from the AFNI python library
# credit to Rick Reynolds, blame to Isaac Schwabacher

import sys

sys.path.append("/apps/x86_64_sci6/AFNI_latest")

import afni_util as UTIL
import lib_afni1D as LAD

class AfniXmat(LAD.Afni1D):
"""Extends Afni1D to allow writing metadata.

lib_afni1D.Afni1D inexplicably disallows writing all the metadata it
collects from the .xmat.1D files it reads. This extension fixes that."""

def write(self, fname, sep=" ", overwrite=0):
"""write the data to a new .1D file

fname : filename is required
sep : separator between elements
overwrite : whether to allow overwrite

return status"""

if self.verb > 2: print '-- Afni1D write to %s, o=%d'%(fname,overwrite)

if not self.ready:
print "** Afni1D not ready for write to '%s'" % fname
return 1
if not fname:
print "** missing filename for write"
return 1

if fname == '-' or fname == 'stdout': fp = sys.stdout
else:
# normal file name: check existence and open
if os.path.exists(fname) and not overwrite:
print "** output file '%s' exists and 'overwrite' not set..."%fname
return 1

fp = open(fname, 'w')
if not fp:
print "** failed to open '%s' for writing" % fname
return 1

fp.write('# <matrix\n')
fp.write('# ni_type = "%d*double"\n' % self.nvec)
fp.write('# ni_dimen = "%d"\n' % self.nt)
fp.write(''.join(['# ColumnLabels = "', ' ; '.join(self.labels), \
'"\n']))
fp.write(''.join(['# ColumnGroups = "', \
UTIL.encode_1D_ints(self.groups), '"\n']))
fp.write('# RowTR = "%f"\n' % self.tr)
fp.write(''.join(['# GoodList = "', \
UTIL.encode_1D_ints(self.goodlist), '"\n']))
fp.write('# NRowFull = "%d"\n' % self.nrowfull)
fp.write(''.join(['# RunStart = "', \
UTIL.encode_1D_ints(self.runstart), '"\n']))
stims = [(self.groups.index(i), \
self.nvec - 1 - list(reversed(self.groups)).index(i)) \
for i in sorted(list(set(self.groups))) if i > 0]
fp.write(''.join(['# Nstim = "%d"\n' % len(stims)]))
fp.write(''.join(['# StimBots = "', UTIL.encode_1D_ints( \
[i[0] for i in stims]), '"\n']))
fp.write(''.join(['# StimTops = "', UTIL.encode_1D_ints( \
[i[1] for i in stims]), '"\n']))
fp.write(''.join(['# StimLabels = "', ' ; '.join([l[:l.rfind('#')] \
for l in [self.labels for (i,j) in stims] ]), '"\n']))
fp.write(''.join(['# CommandLine = "%s ; ' % self.command, \
' '.join(sys.argv), '"\n']))
fp.write('# >\n')

for row in range(self.nt):
fp.write(sep.join(['%g' % self.mat[row] for i in range(self.nvec)]))
fp.write('\n')

fp.write('# </matrix>\n')

fp.close()

return 0

xmat = AfniXmat(sys.argv[1])
nonzerocols = [i for i in xrange(len(xmat.mat)) if min(xmat.mat) < -.00005 \
or max(xmat.mat) > .00005] # epsilon = .00005 in 3dDeconvolve
xmat.reduce_by_vec_list(nonzerocols)
xmat.write('-')
Subject Author Posted

Auto-dropping all zero regressors

Andrew Brooks July 10, 2012 02:31PM

Re: Auto-dropping all zero regressors

Andrew Brooks July 10, 2012 02:53PM

Re: Auto-dropping all zero regressors

Isaac Schwabacher July 10, 2012 03:11PM

Re: Auto-dropping all zero regressors

Andrew Brooks July 10, 2012 04:14PM

Re: Auto-dropping all zero regressors

Isaac Schwabacher July 10, 2012 04:29PM

Re: Auto-dropping all zero regressors

Andrew Brooks July 10, 2012 05:57PM