[MINC-users] Calculating an FDR threshold

Jason Lerch jason at phenogenomics.ca
Wed Sep 28 07:45:26 EDT 2011


There are routines inside RMINC (a library in R) that can calculate the FDR for MINC volumes - mincFDR.default in particular.

If you have R but not RMINC (and can't be bothered to install it), here's a moderately quick workaround (all this assumes a two tailed t-test; it's easy to adjust for another type of test):

mincextract -ascii statsvol.mnc > stasvol.txt
R
> vol <- read.table('statsvol.txt')$V1
> pvals <- 2* pt(-abs(vol), df) # where df is your degrees of freedom
> qvals <- p.adjust(pvals, 'fdr')
> max(pvals[qvals<0.05]) # the p value corresponding to a 5% FDR
> qt(max(pvals[qvals<0.05])/2, df) # the t value corresponding to a 5% FDR

Using a mask is left as an exercise to the reader ;-)

Jason

On 2011-09-27, at 11:23 PM, Mishkin Derakhshan <mishkind at gmail.com> wrote:

> Hi Jennika,
> 
> I've never used glim_image or fMRIstat, but I have computed FDR
> thresholds using R and the mni.cortical.statistics library.
> 
> If you can get the VBM t-stat map that is generated from glim_image,
> into a simple text file listing each t-stat per row, then you can use
> the following:
> 
> 1. Install R and library
> on ubuntu: sudo apt-get install r-base-core
> then download the mni.cortical.statistics library to /tmp/
> http://www.bic.mni.mcgill.ca/users/mishkin/files/mni.cortical.statistics.0.9.0_patched.tar.gz
> cd /tmp
> tar xvzf mni.cortical.statistics.0.9.0_patched.tar.gz
> mkdir /home/mishkin/R_LIBS (or wherever you want to keep your R libraries)
> R CMD INSTALL --library=/home/mishkin/R_LIBS/ mni.cortical.statistics
> export R_LIBS=/home/mishkin/R_LIBS
> 
> 2. Calculate the FDR in R
> R
> library(mni.cortical.statistics)
> tstats =scan("my_file_of_tstats.txt")
> q <- mni.compute.FDR(tstats, df=12, fdr=0.05)
> write(q,"fdr_threshold_output.txt",1)
> 
> where:
> "my_file_of_tstats.txt" looks like this:
> -2.23
> -2.24
> ...
> -1.93
> 0.32
> 
> "df = 12" means 12 degrees of freedom. usually the N of each group -2
> "fdr = 0.05" means use an fdr threshold of 0.05
> 
> "fdr_threshold_output.txt" will be an output file containing a number
> that indicates the fdr threshold t-stat:
> -2.1
> 
> So, this would mean that all t-stats less than -2.1 would be deemed
> significant using an FDR of 0.05.
> 
> hth,
> mishkin
> 
> 
> On Tue, Sep 27, 2011 at 9:04 PM, Jennika Soles <jennika.soles at gmail.com> wrote:
>> Hi All,
>> 
>> Does anyone know if there is a particular method of calculating an FDR
>> threshold for a VBM t-stat map generated using glim_image? Would using the
>> tool in fMRIstat be acceptable?
>> 
>> Thanks!
>> _______________________________________________
>> MINC-users at bic.mni.mcgill.ca
>> http://www.bic.mni.mcgill.ca/mailman/listinfo/minc-users
>> 
> _______________________________________________
> MINC-users at bic.mni.mcgill.ca
> http://www.bic.mni.mcgill.ca/mailman/listinfo/minc-users


More information about the MINC-users mailing list