After further investigation this morning, we have found that the problem with 3dDeconvolve occurs in gcc with any level of optimization in the eis_svd.c function, but apparently only on Intel-like Linux systems. I have recompiled all 6 binary packages now with all optimization turned off for this function, and these binaries (and the source) are available for download as of now (binaries you fetched this morning may not be good, sorry). When you run 3dDeconvolve, the latest revision date should be today (19 Aug 2004). Anything after 19 Jul 2004 and before today should be considered suspicious.
The problem does not occur with all matrices (which are assembled from the
-stim_files), but I have no way to predict when the problem does occur. As a diagnostic for this problem, 3dDeconvolve has a new output line, which appears right after the condition number:
++ (438x96) Matrix condition [XtX]: 60.1337
++ Matrix inverse average error = 8.0531e-17
This "inverse average error" is the average absolute difference between the product of the
X matrix pseudoinverse and the
X matrix itself, minus the
I matrix. Ideally, this value should be zero. Roundoff error in the calculations, in double precision, is about 1 part in 10^16, as in the example above (in single precision -- 3dDeconvolve_f -- should be about 1 part in 10^7).
You can test your runs from the last 3-4 weeks using the
-nodata option, setting things up the same otherwise, and looking at the "inverse average error" (e.g., with grep). Any value significantly larger than roundoff is bad bad bad. In the case that started this imbroglio yesterday afternoon, the "inverse average error" is about 0.01.
This problem may occur with regression or with deconvolution -- the same C functions are used to setup and solve the equations for both cases.
The problem did not occur on the Mac OS X (gcc 3.3) or Solaris (Sun's C) platforms, but did occur with gcc 2.96, gcc 3.2, and gcc 3.3 on Intel-ish platforms. It did not occur on the Opteron compiled with gcc 3.3 in 64 bit mode. The problem is
not related to the Legendre polynomials.
Again, I'm sorry for this, esp. for the extra work you have to do to check the last few weeks runs --
bob cox
P.S. The 6 binary packages we keep up-to-date are
* linux_gcc32 = compiled with gcc 3.2.2 on a RedHat 9 system
* linux_glibc22 = compiled with gcc 2.96 on a RedHat 7.3 system
* macosx_10.3_G5 = compiled with gcc 3.3 on a Mac OS X 10.3 G5
* macosx_10.3 = compiled with gcc 3.3 on a Mac OS X 10.3 G4
* sgi10k_6.5_gcc = compiled with gcc 3.01 on an SGI
* solaris29_suncc = compiled with SUNWspro on a Solaris 2.9 system
The NIH machines on which I remotely install, and to which I still have access (some of you no longer let me login) have been updated.