Millisecond Forums

gaussian ISI

https://forums.millisecond.com/Topic8431.aspx

By Mela - 7/17/2012

Hi all,


lets say we want an experiment where the interval between stimuli is variable - one simple way to achieve this is to use the rand function and set the boundaries (e.g., rand(200, 1000). But what if the requirement is for the ISI to be randomly selected from a gaussian distribution with known mean & SD - any possibilities here?


I know the monkey allows such style of distributions, but using the normal(x, y) doesn't appear to work alongside the post/pre trial commands...


Any solutions/ideas?


Thanks,


Mela


By Dave - 7/17/2012

You can implement this yourself via Inquisit's expressions syntax. See e.g. https://www.millisecond.com/forums/Topic5419.aspx#5419 for somewhat related* code (standard normal CDF and its inverse). In short, pick a suitable algorithm to generate random Gaussian variates from the literature (there are many) and express it in Inquisit's terms.


Regards,


~Dave


* e.g. https://en.wikipedia.org/wiki/Inverse_transform_sampling

By Mela - 7/17/2012

Super, thanks Dave. Certainly a viable solution!


Would be nice to have a simple command, of course (can't be hard to implement given the element is there for the monkey!).

By Dave - 7/17/2012

Would be nice to have a simple command


It's on the list, but not implemented yet as there are other features with higher priority.

By Mela - 7/17/2012

Sure. Good to hear it's already on the radar, though!


Thanks for the help, Dave (you're a star!).

By Mela - 7/18/2012

OK, so had  a little think about the nice code and info you linked to and here's where I think I've got to...


so the idea would be to...


1. have a random selection from values ranging between 0 and 1 (which is a uniform distribution) - this gives the u value (as defined on wiki). Essentially the probability value.


2. use the z_p expression you linked to (or similar) which uses the above u value and calculates its associated z score


3. Use z score, along with defined distribution mean and SD to calculate ISI


So if this is implemented ontrialbegin it would produce a variable ISI which is effectively a random value from the defined distribution. That's a pretty elegant solution. I haven't implemented it yet, but certainly something I should get to work (given the above is sound).


By Dave - 7/18/2012


1. have a random selection from values ranging between 0 and 1 (which is a uniform distribution) - this gives the u value (as defined on wiki). Essentially the probability value.


2. use the z_p expression you linked to (or similar) which uses the above u value and calculates its associated z score


3. Use z score, along with defined distribution mean and SD to calculate ISI



Yep, that's the way to do it.


EDIT: For the sake of posterity, here's syntax to generate samples from a normal distribution with specified mean mu and standard deviation sigma:



/* values.nsamples:    number of random variates (samples) to draw from distribution
/* values.mu:          desired (theoretical) mean of normal distribution
/* values.sigma:       desired (theoretical) standard deviation of normal distribution
/* expressions.mx:     observed (empirical) mean of random variates
/* expressions.sdx:    observed (empirical) standard deviation of random variates

<expressions>
/ z = if(values.u>0.5)(sqrt(-ln(1-pow(2*values.u-1,2))/sqrt(m_pi/8))) else
    (-sqrt(-ln(1-pow(2*values.u-1,2))/sqrt(m_pi/8)))
/ x = (values.sigma*expressions.z)+values.mu
/ mx = (values.sumx/values.nx)
/ sdx = sqrt((values.ssumx-(values.nx *(expressions.mx*expressions.mx)))/(values.nx-1))
</expressions>

<values>
/ nsamples = 2000
/ mu = 1500
/ sigma = 500
/ u = 0
/ sumx = 0
/ ssumx = 0
/ nx = 0
</values>

<text info>
/ items = ("Input: mu = <%values.mu%> | sigma = <%values.sigma%> | n = <%values.nsamples%>
Output: m(x) = <%expressions.mx%>| sd(x) = <%expressions.sdx%> | n(x) = <%values.nx%>")
/ size = (90%,90%)
/ erase = false
/ vjustify = center
</text>

<trial sample>
/ ontrialbegin = [values.u=rand(0,1); values.sumx=values.sumx+expressions.x; values.ssumx=values.ssumx+(expressions.x*expressions.x); values.nx=values.nx+1]
/ validresponse = (noresponse)
/ trialduration = 0
/ branch = [if(block.myblock.trialcount<values.nsamples)trial.sample else trial.summary]
</trial>

<trial summary>
/ stimulusframes = [1=info]
/ validresponse = (57)
/ recorddata = false
</trial>

<block myblock>
/ trials = [1=sample]
</block>

<data>
/ columns = [trialnum,expressions.x,expressions.mx,expressions.sdx]
/ separatefiles = true
</data>

<defaults>
/ windowsize = (640px,480px)
/ fontstyle = ("Verdana", 2%, true)
</defaults>


To check the output and verify that it's indeed generating the desired normal samples, generate a histogram of expressions.x, examine its Q-Q plot and/or submit to a statistical test for normality (Kolmogorov-Smirnov, Shapiro-Wilk, etc.).


Regards,


~Dave

By Mela - 7/18/2012

Cheers!


May Cthulu eat your house last :)

By Dave - 7/18/2012