6

I want to be able to compute following integral using numpys trapz function

numpy.trapz([-1, 1]) # returns 0

But I don't want to allow negative areas. Is there an efficient way to do this or do I have to look for the minimum point and transform the array by hand?

Does numpy.trapz(numpy.abs([-1, 1])) make sense?

Saullo G. P. Castro
  • 56,802
  • 26
  • 179
  • 234
anopheles
  • 434
  • 1
  • 9
  • 19
  • yes, `np.trapz(np.abs([-1, 1]))` makes sense, but the integral value will be, in this case, the double the value of the positive area, so you have to divide it by the proper ratio, which in this case is `2`... – Saullo G. P. Castro Oct 28 '13 at 16:45
  • When you say "negative areas" do you mean `y` values less then zero or negative contributions to the integration volume? – Daniel Oct 28 '13 at 17:17
  • I want to have a the area of the graph which is independent of the coordinate system which it is in. I know that this is not possible with non closed surfaces. I think taking the abs does the trick. @SaulloCastro How do you compensate the ratio in the general case? – anopheles Oct 28 '13 at 17:55
  • 1
    @anopheles - The area _under_ a curve is never independent of the coordinate system. Are you asking for the area of a specific polygon? If so, you don't want `trapz`. – Joe Kington Oct 28 '13 at 18:29
  • @JoeKington Please read my comment above yours. I'm aware of that. – anopheles Oct 28 '13 at 19:20
  • 1
    @anopheles - In that case, you don't want 1D integration, which is way `trapz` does. It sounds like you want to calculate the area of a polygon in 2D space. It's only a few lines of code (e.g. http://stackoverflow.com/questions/451426/how-do-i-calculate-the-surface-area-of-a-2d-polygon), but if you're doing this sort of thing a lot, have a look at `shapely`. It's a very nice library for geometric operations (e.g. area would just be `Polygon([coords]).area`). – Joe Kington Oct 28 '13 at 19:31

1 Answers1

3

If you want to discard negative contributions to the integrated area we can simply grab the np.trapz source code and rewrite it:

def abstrapz(y, x=None, dx=1.0):
    y = np.asanyarray(y)
    if x is None:
        d = dx
    else:
        x = np.asanyarray(x)
        d = np.diff(x)
    ret = (d * (y[1:] +y[:-1]) / 2.0)
    return ret[ret>0].sum()  #The important line

A quick test:

np.trapz([-1,0,1])
0.0

abstrapz([-1,0,1])
0.5

If you just want to avoid areas where y is less then zero simply mask 'y' values less then zero to zero:

arr = np.array([-2,-1,0.5,1,2,1,0,5,3,0])
np.trapz(arr)
10.5

arr[arr<0] = 0
np.trapz(arr)
12.5

This isn't the best way to do it, but it is a good approximation. If this is what you mean I can update this.

I had to change your example slightly as trapz([-1,1]) will always return 0 by definition. We do remove some functionality this way, if you need to do this on multidimensional arrays it is easy to add it back in.

Daniel
  • 19,179
  • 7
  • 60
  • 74
  • By doing ret[ret>0] don't you just discard every negative contribution? I want to have an area in the traditional sense, i.e. treat every subarea positive. In your example, an array with only negative values has a surface of zero. This is not what I want. – anopheles Oct 28 '13 at 17:45