Re: [ferret_users] event statistics -- length of contiguous "spell"

[Steve Hankin <steven.c.hankin@xxxxxxxx>]

Hi Jagadish,

Here is an example of using the "closest index above" and "closest index below" filters in Ferret to compute the length of contiguous "spells"

define axis/t=1-jan-1900:1-jan-2000:30/units=days t30
let spells = IF RANDU(L[gt=t30]) GT 0.5 THEN 1
let inverseSpells = IF MISSING(spells,2) EQ 2 THEN 1
let closestAbove = inverseSpells[l=@cia]
let closestBelow = inverseSpells[l=@cib]
let rawSize = closestAbove+closestBelow
let contigSize = IF rawSize NE 0 THEN rawSize-1
 
list/l=1:50 spells, inverseSpells, closestAbove, closestBelow, contigSize
             TIME: 17-DEC-1899 00:00 to 26-JAN-1904 00:00
 Column  1: SPELLS is IF RANDU(L[GT=T30]) GT 0.5 THEN 1
 Column  2: INVERSESPELLS is IF MISSING(SPELLS,2) EQ 2 THEN 1
 Column  3: CLOSESTABOVE is INVERSESPELLS[L=@CIA]
 Column  4: CLOSESTBELOW is INVERSESPELLS[L=@CIB]
 Column  5: CONTIGSIZE is IF RAWSIZE NE 0 THEN RAWSIZE-1
                 SPELLS  INVERSE CLOSEST CLOSEST CONTIGSIZE
01-JAN-1900 /  1:   ....   1.000   0.000   0.000    ....
31-JAN-1900 /  2:  1.000    ....   1.000   1.000   1.000
02-MAR-1900 /  3:   ....   1.000   0.000   0.000    ....
01-APR-1900 /  4:  1.000    ....   1.000   1.000   1.000
01-MAY-1900 /  5:   ....   1.000   0.000   0.000    ....
31-MAY-1900 /  6:  1.000    ....   1.000   1.000   1.000
30-JUN-1900 /  7:   ....   1.000   0.000   0.000    ....
30-JUL-1900 /  8:  1.000    ....   1.000   1.000   1.000
29-AUG-1900 /  9:   ....   1.000   0.000   0.000    ....
28-SEP-1900 / 10:   ....   1.000   0.000   0.000    ....
28-OCT-1900 / 11:   ....   1.000   0.000   0.000    ....
27-NOV-1900 / 12:   ....   1.000   0.000   0.000    ....
27-DEC-1900 / 13:  1.000    ....   1.000   1.000   1.000
26-JAN-1901 / 14:   ....   1.000   0.000   0.000    ....
25-FEB-1901 / 15:  1.000    ....   2.000   1.000   2.000
27-MAR-1901 / 16:  1.000    ....   1.000   2.000   2.000
26-APR-1901 / 17:   ....   1.000   0.000   0.000    ....
26-MAY-1901 / 18:   ....   1.000   0.000   0.000    ....
25-JUN-1901 / 19:  1.000    ....   2.000   1.000   2.000
25-JUL-1901 / 20:  1.000    ....   1.000   2.000   2.000
24-AUG-1901 / 21:   ....   1.000   0.000   0.000    ....
23-SEP-1901 / 22:   ....   1.000   0.000   0.000    ....
23-OCT-1901 / 23:   ....   1.000   0.000   0.000    ....
22-NOV-1901 / 24:   ....   1.000   0.000   0.000    ....
22-DEC-1901 / 25:   ....   1.000   0.000   0.000    ....
21-JAN-1902 / 26:  1.000    ....   2.000   1.000   2.000
20-FEB-1902 / 27:  1.000    ....   1.000   2.000   2.000
22-MAR-1902 / 28:   ....   1.000   0.000   0.000    ....
21-APR-1902 / 29:  1.000    ....   1.000   1.000   1.000
21-MAY-1902 / 30:   ....   1.000   0.000   0.000    ....
20-JUN-1902 / 31:  1.000    ....   3.000   1.000   3.000
20-JUL-1902 / 32:  1.000    ....   2.000   2.000   3.000
19-AUG-1902 / 33:  1.000    ....   1.000   3.000   3.000
18-SEP-1902 / 34:   ....   1.000   0.000   0.000    ....
18-OCT-1902 / 35:   ....   1.000   0.000   0.000    ....
17-NOV-1902 / 36:  1.000    ....   2.000   1.000   2.000
17-DEC-1902 / 37:  1.000    ....   1.000   2.000   2.000
16-JAN-1903 / 38:   ....   1.000   0.000   0.000    ....
15-FEB-1903 / 39:   ....   1.000   0.000   0.000    ....

On 5/14/2014 8:56 AM, Steve Hankin wrote:

On 5/13/2014 5:19 PM, 'jagadish karmacharya' via _OAR PMEL Ferret Users wrote:

Dear ferret users,

I have a 1 dimensional time series with intermittent values equal to 1 and missing values in between (1's represent spells of certain events). I want to calculate number of spells of of 1's and other stats such as mean and medium duration of spells of 1's.

Hi Jagadish,

Following up on the ideas that Hein has begun ...

Your description is not sufficiently specific on the definition of a "spell" to give a precise answer on how (or whether) the tools available in Ferret can do what you are after.

• if the definition of a spell is that it a threshold on the *density* of 1's found per unit length of the time series, then the @SBX (running average) filter will give you something very close to what you are after
• e.g. yes? let tdensity = mytimeseries[l=@sbx:11]   
         yes? PLOT tdensity    ! have a look ...
         yes? let myspells = IF tdensity gt threshold THEN 1
  
• if the definition of a spell is that it must be a strictly contiguous sequence of 1's, it will take more imagination to see how to apply the available Ferret tools.    Here's one approach to identifying the contiguous spells using existing tools:
1. use the tdensity approach above varying the box width to identify the longest spell. 
   yes? let tdensity = mytimeseries[l=@sbx:11]
   yes? let elevenOrLonger = IF tdensity EQ 1 THEN 1
   yes? LIST elevenOrLonger[L=@SUM]     ! if answer is zero, there are no spells this long
   
2. define a new time series in which you mask out the blocks of contiguous 1s just found
   Iterate the procedure to find number and location of spells of shorter and shorter duration.  (Note only odd-number length spells can be found this way.)
   
• Or ... as Hein says, you can write raw code to perform any custom calculation you want by using either Python in pyFerret or defining your own external function in vanilla Ferret.  ==> Then contribute the idea back to the Ferret user's list!  A filter that, at each point, measured the length of contiguous block of non-missing values surrounding it, would be a nice, generic Ferret tool.  (Should we add such a filter into Ferret?  It would be a fairly easy addition:   e.g. LET contiguous_spell_length = mytimeseries[l=@cntg])
• I cant help wondering if a clever use of the nearest index above and nearest index below filters might be used to build a contiguous length filter using a few LET definitions ...
  

    - Steve

I guess @evnt transformation would be helpful for this, in particular defining ID number as outlined here http://www.pmel.noaa.gov/maillists/tmap/ferret_users/fu_2013/msg00277.html would provide an solution. I made an attempt accordingly but haven't been able to solve the problem.

I would appreciate any lead on this.

Thanks,

Jagadish