MathGroup Archive 2008

[Date Index] [Thread Index] [Author Index]

Search the Archive

Re: Re: Problem with BinCounts

  • To: mathgroup at smc.vnet.net
  • Subject: [mg90950] Re: [mg90927] Re: Problem with BinCounts
  • From: Darren Glosemeyer <darreng at wolfram.com>
  • Date: Thu, 31 Jul 2008 02:57:16 -0400 (EDT)
  • References: <g6mas3$i0d$1@smc.vnet.net> <200807300751.DAA17667@smc.vnet.net>

Valeri Astanoff wrote:
> On 29 juil, 07:46, Mark Teagarden <Mark.Teagar... at UTSA.EDU> wrote:
>   
>> Hi,
>>
>> I'm trying to do an autocorrelation analysis on an artificially generated
>> neuronal spike train, but BinCounts is behaving strangely.
>>
>> I first generate a list of spike times for an oscillator with a mean
>> frequency of 1 Hz, plus some gaussian jitter:
>>
>> st = NestList[(1 + RandomReal[NormalDistribution[0, 0.1]] + # &), 0,100=
>>     
> 0];
>   
>> Then I do the autocorrelation like this:
>>
>> BinCounts[st, {# - 2.00, # + 2.00, 0.01}] & /@ st;
>>
>> ac = Last[Accumulate[%]]/Length[%];
>>
>> The output of BinCounts ought to be a list of length st, with each list
>> element having a fixed length of 400.  And yet I get error messages lik=
>>     
> e so:
>   
>> Thread::tdlen: Objects of unequal length in
>> \
>> {13,11,15,9,12,13,14,13,9,6,<<390>>}+{0,1,0,0,0,0,0,0,0,0,<<389>>}
>> \
>> cannot be combined. >>
>>
>> And sure enough, if I run Dimensions/@ the output of the BinCounts line, =
>>     
> I
>   
>> get one or two elements with a length of 399.  This happens completely =
>>     
> at
>   
>> random, and I don't understand why.  It seems to me that BinCounts ough=
>>     
> t to
>   
>> give me a fixed output length, regardless of what the results of the
>> operation actually were.  These shorter elements do not occur near the
>> beginning or end, either, so it's not an edge effect.  Furthermore, the
>> values of st corresponding to the shorter elements of ac don't have any
>> unusual properties like being whole numbers or anything.  I will add th=
>>     
> at
>   
>> sometimes the code works correctly, and it's not obvious why.
>>
>> Has anyone else encountered such behavior before?  What did you do abou=
>>     
> t it?
>   
>
> Good day,
>
> All I can say is that if you rationalize 'st'
> and all parameters then it doesn't occur :
>
> st=NestList[(1+RandomReal[NormalDistribution[0,0.1]]+#&),
> 0,1000]//Rationalize[#,10^-6]&;
>
> BinCounts[st, {# - 2, # + 2, 1/100}] & /@ st;
>
>
> V.Astanoff
>
>   

It's clear that this is a consequence of inexact arithmetic, but not yet 
clear if it is avoidable in general. Basically, a small amount of 
numerical error is enough to cause a bin to not be included (effectively 
something like a step going slightly beyond the end point in the inexact 
representation, though not with exact numbers) for some starting values.

A way around this is to shift the data being binned rather than the bins 
themselves, as we see here in 100 trials

In[1]:= Table[st =
           NestList[(1 + RandomReal[NormalDistribution[0, 0.1]] + # &), 0,
            1000];
          bc2 = BinCounts[st - #, {-2.00, 2.00, 0.01}] & /@ st;
          Map[Length, bc2] // Union, {100}] // Union

Out[1]= {{400}}


and better still to do this and use an exact binning spec to avoid any 
numerical error in the bins.


In[2]:= Table[st =
           NestList[(1 + RandomReal[NormalDistribution[0, 0.1]] + # &), 0,
            1000];
          bc2 = BinCounts[st - #, {-2, 2, 1/100}] & /@ st;
          Map[Length, bc2] // Union, {100}] // Union

Out[2]= {{400}}


Darren Glosemeyer
Wolfram Research


  • Prev by Date: Re: Re: bug? f'[x]'
  • Next by Date: rotate axes labels with Plot3D?
  • Previous by thread: Re: Problem with BinCounts
  • Next by thread: remove higher orders terms from mathemtica's result