Calculating mode using rangestat command

[Ruzica Savcic]

Is there a way to find the mode of a variable within a specified window of time using the rangestat command?
I've never written a function in Mata before and after reading a few threads on this topic, I came up with this

mata:
     real rowvector mode(real matrix X) {
        return(mode(X), max)
    }
end

rangestat (mode) my_var, interval(time -10 10)

The programme returns error 3499, "function not found". How do you specify mode in Mata? I tried searching in Mata manual, but the word 'mode' is used in a different sense there.

[Robert Picard]

I don't think there's a mode() function in Mata. The options are to write one and use it with rangestat or use rangerun and use Stata's egen mode(). The latter is much simpler to do:

clear all
set obs 30
gen time = _n
set seed 312
gen my_var = runiformint(1,20)

program my_mode
    egen wanted = mode(my_var), min
end
rangerun my_mode, interval(time -2 2)

[Nick Cox]

There is no built-in for mode in rangestat (SSC) and no such Mata function bundled with Stata.

That's not fatal, but you'd need to write one.

Let's answer a question with a question: what's the mode of 1,1,2,2,3,3? If you can specify your rules, you're half-way to an answer.

It would be possible to do it in other ways. Here's one. I take hsmode from SSC and use the Mata code. (The author is me, so no permission or originality issues.)

webuse grunfeld , clear 

mata

mata clear 

real hsmode(real colvector y) 
{ 
        real scalar n, mode 

        n = rows(y) 

        if (n <= 3) { 
                mode = report(y)
        }
        else { 
                real matrix diff 
                real scalar n2, i1, i2, nties

                i1 = 1 
                _sort(y, 1) 

                while (n > 3) { 
                        n2 = floor(n/2) 
                        n1 = n - n2
                        i2 = i1 + n2 
                        diff = y[i2 .. (i2 + n1 - 1)] - y[i1 .. (i1 + n1 - 1)] 
                        diff = diff, (i1 .. (i1 + n1 - 1))' 
                        _sort(diff, (1,2))
                        nties = colsum(diff[,1] :== diff[1,1]) 
                        i1 = diff[ceil(nties / 2), 2] 
                        n = n1 
                }

                mode = report(y[i1 .. (i1 + n1 - 1)]) 
        }

        return(mode) 
}       

real report(real colvector y) 
{ 
        if (rows(y) == 1) { 
                return(y[1])
        }
        else if (rows(y) == 2) {
                return((y[1] + y[2]) / 2) 
        }
        else if (rows(y) == 3) { 
                _sort(y, 1) 
                if ((y[2] - y[1]) < (y[3] - y[2])) { 
                        return((y[1] + y[2]) / 2)
                }
                else if ((y[2] - y[1]) > (y[3] - y[2])) { 
                        return((y[3] + y[2]) / 2)
                }
                else return(y[2]) 
        } 
}       

end 

rangestat (hsmode) invest, int(year -5 5) by(company) 

list company year invest invest_count hsmode1 in 1/20, sepby(invest_count)

     +------------------------------------------------+
     | company   year   invest   invest~t     hsmode1 |
     |------------------------------------------------|
  1. |       1   1935    317.6          6       324.2 |
     |------------------------------------------------|
  2. |       1   1936    391.8          7       324.2 |
     |------------------------------------------------|
  3. |       1   1937    410.6          8       429.3 |
     |------------------------------------------------|
  4. |       1   1938    257.7          9   454.60001 |
     |------------------------------------------------|
  5. |       1   1939    330.8         10   454.60001 |
     |------------------------------------------------|
  6. |       1   1940    461.2         11       505.8 |
  7. |       1   1941      512         11       505.8 |
  8. |       1   1942      448         11   554.35001 |
  9. |       1   1943    499.6         11   554.35001 |
 10. |       1   1944    547.5         11   558.14999 |
 11. |       1   1945    561.2         11   558.14999 |
 12. |       1   1946    688.1         11   558.14999 |
 13. |       1   1947    568.9         11   558.14999 |
 14. |       1   1948    529.2         11   558.14999 |
 15. |       1   1949    555.1         11   558.14999 |
     |------------------------------------------------|
 16. |       1   1950    642.9         10   558.14999 |
     |------------------------------------------------|
 17. |       1   1951    755.9          9         562 |
     |------------------------------------------------|
 18. |       1   1952    891.2          8   542.14999 |
     |------------------------------------------------|
 19. |       1   1953   1304.4          7   542.14999 |
     |------------------------------------------------|
 20. |       1   1954   1486.7          6         599 |
     +------------------------------------------------+

.

[Ruzica Savcic]

Thank you, Robert and Nick!
rangerun, of course!
My exercise is such that in case of polymodal variables, I would opt for the highest value, hence

 
 egen wanted = mode(my_var), max

[Harry Kaj]

I ran the code below as suggested, but I was wondering if there is a way to include how frequent the mode was; whether the mode occurred 5 times or 7 times or 10 times etc?

I don't think there's a mode() function in Mata. The options are to write one and use it with rangestat or use rangerun and use Stata's egen mode(). The latter is much simpler to do:

clear all
set obs 30
gen time = _n
set seed 312
gen my_var = runiformint(1,20)

program my_mode
egen wanted = mode(my_var), min
end
rangerun my_mode, interval(time -2 2)

[Nick Cox]

Naturally what to do is to extend the program written for rangerun

clear all
set obs 30
gen time = _n
set seed 312
gen my_var = runiformint(1,20)

program my_mode
egen wanted = mode(my_var), min
count if wanted == my_var 
gen wanted2 = r(N)
end

rangerun my_mode, interval(time -2 2)

[Harry Kaj]

Of course! thank you, Nick