Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Re-implement stat_surrogate_pvals #299

Merged
merged 3 commits into from
May 17, 2021
Merged

Re-implement stat_surrogate_pvals #299

merged 3 commits into from
May 17, 2021

Conversation

tstenner
Copy link
Contributor

Taken from sccn/SIFT#2:

stat_surrogateStats forwards the statTest.tail parameter to stat_surrogate_pvals which does a two-tailed test when statTest.tail is both and a right-tailed test for all other inputs.

Despite this, stat_surrogateStats prints

'This is a left-sided test for significant differences between conditions'

even though the p values are for a right sided test.

This PR reimplements the stat_surrogate_pvals function.
The main differences are:

  • the algorithm is O(n) instead of O(log(n)) (benchmark at the bottom)
  • it's a bit more precise, i.e. it calculates p=sum(t ≥ T)/n instead of p=(sum(t ≥ T)+.5)/(n+1)
  • 'left' is accepted as tail parameter
  • an error is raised when tail is an unexpected value instead of assuming a right-sided test
>> rng('default');
>> distribution = rand(5,7,200);
>> observed = reshape(linspace(0, 1, 35), 5, 7);
>> % new implemantation
>> stat_surrogate_pvals(distribution, observed, 'right')
ans =
    1.0000    0.8350    0.7350    0.5500    0.4650    0.2900    0.1350
    0.9750    0.8150    0.7550    0.5250    0.3250    0.2850    0.0650
    0.9550    0.7750    0.6750    0.5100    0.3600    0.2150    0.0600
    0.8750    0.7900    0.6050    0.4350    0.3400    0.1800    0.0200
    0.8650    0.7200    0.5800    0.4200    0.3100    0.1100         0
>> % old implementation
>> stat_surrogate_pvals_(distribution, observed, 'right')
ans =
    0.9975    0.8333    0.7338    0.5498    0.4652    0.2910    0.1368
    0.9726    0.8134    0.7537    0.5249    0.3259    0.2861    0.0672
    0.9527    0.7736    0.6741    0.5100    0.3607    0.2164    0.0622
    0.8731    0.7886    0.6045    0.4353    0.3408    0.1816    0.0224
    0.8632    0.7189    0.5796    0.4204    0.3109    0.1119    0.0025
>> stat_surrogate_pvals(distribution, observed, 'both')
ans =
         0    0.3300    0.5300    0.9000    0.9300    0.5800    0.2700
    0.0500    0.3700    0.4900    0.9500    0.6500    0.5700    0.1300
    0.0900    0.4500    0.6500    0.9800    0.7200    0.4300    0.1200
    0.2500    0.4200    0.7900    0.8700    0.6800    0.3600    0.0400
    0.2700    0.5600    0.8400    0.8400    0.6200    0.2200         0
>> stat_surrogate_pvals_(distribution, observed, 'both')
ans =
    0.0050    0.3333    0.5323    0.9005    0.9303    0.5821    0.2736
    0.0547    0.3731    0.4925    0.9502    0.6517    0.5721    0.1343
    0.0945    0.4527    0.6517    0.9801    0.7214    0.4328    0.1244
    0.2537    0.4229    0.7910    0.8706    0.6816    0.3632    0.0448
    0.2736    0.5622    0.8408    0.8408    0.6219    0.2239    0.0050
>> % are row / column vectors handled properly?
>> stat_surrogate_pvals(1:100, 4.5, 'left')
ans =
    0.0400
>> stat_surrogate_pvals(reshape(1:100, [], 1), 4.5, 'left')
ans =
    0.0400

Speed:

>> timeit(@() stat_surrogate_pvals(distribution, observed, 'right'))
ans =
   1.7884e-05
>> timeit(@() stat_surrogate_pvals_(distribution, observed, 'right'))
ans =
   1.2996e-04
>> distribution = rand(500,700,200);
>> observed = reshape(linspace(0,1, 350000), 500,700);
>> timeit(@() stat_surrogate_pvals(distribution, observed, 'right'))
ans =
    0.0485
>> timeit(@() stat_surrogate_pvals_(distribution, observed, 'right'))
ans =
    1.1268

@arnodelorme arnodelorme merged commit 7b1d955 into sccn:develop May 17, 2021
@tstenner tstenner deleted the surrogate_pvals branch May 17, 2021 17:23
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

2 participants