How can I take the cumulative sum of a 2D distribution, starting from the highest value in the distribution?

Background: my 2D distribution is the probability of different combinations of wave heights and wave periods based on measurements made at an ocean site. Assume the distribution shows that the most common wave (highest probability, 1%) has a height of 3 meters and a period of 8 seconds. This combination represents the "peak" of a "mountain" in the distribution, i.e. adjacent combinations of wave height/period occur often but less frequently (say 0.1%), and combinations further "downhill" occur even less frequently (say 0.01%). Non-zero probability is spread out over many bins of wave height and wave period so it still totals to 1.

I want to be able to identify the sea states that occur 10% of the time, i.e. make a contour that encloses the height/period combinations that together add up to a probability of 10%. In other words, I can look at the contour and say "10% of the time, wave heights are between 1 and 4 meters and periods are between 6 and 10 seconds". I believe I have to contour the cumulative sum of the 2D distribution starting at the "peak" and working "downhill". What is confusing is that downhill extends in 2 directions since I am summing over 2 variables, height and period.

Any advice is helpful since I am having trouble wrapping my mind around this.

  • $\begingroup$ It sounds like you're after the coverage of a highest density region. $\endgroup$ – Glen_b Aug 13 '16 at 5:55

Here is one approach I came up with, but if there's a more elegant/accurate approach please let me know.

  1. Re-sample the 2D histogram to very small bins of wave height and period, and normalize so the sum of the resampled probability is 1. Smaller bins make the summation below more accurate.
  2. Loop through different probability values. For each probability value, plot the contour of that probability and identify enclosed bins of wave height and period (I am doing this in Matlab and used the approach in this answer).
  3. Sum the probability in the enclosed bins and store the value.
  4. After the loop is finished, interpolate to find the contour that would enclose 10%. If x=[summed probability from step 3] and y=[contoured probability value from step 2], interpolate on xi=[desired summed probability], here 10%, to find yi=[contour enclosing desired summed probability].
  5. Delete all prior tested contours from the plot, then plot the final contour of yi which encloses the height/period combinations that happen 10% of the time.
| cite | improve this answer | |
  • 4
    $\begingroup$ This is the right approach, but the solution can be obtained in a much easier fashion: having rasterized a density estimate, all you need to do is invert the cumulative sum of the values on the raster. This is much faster and simpler than looping over values and plotting contours. (This question was asked several years ago on the GIS site; answers include working R code, which requires all of six lines to find and display the solution.) A closely related question appears at stats.stackexchange.com/questions/63447. $\endgroup$ – whuber Aug 18 '14 at 14:25

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.