Search code examples
arraysmatlabmultidimensional-arrayplotprobability-density

Obtain P(x)=integrate W(v,x)dv, from the matrix A(Vvalues,Xvalues,W(v,x)) in Matlab?


OK, first to clarify: I am asking for an algorithm that I could use for my problem please (however to get a Matlab code in addition would be nice). I didn't post any code here, because I don't know how to write an algorithm for this.

The MATLAB matrix A has three columns, v,x,W(v,x)), for v in the range

v_min <= v <= v_max

and x in the range

x_min <= x <= x_max 

The first two columns have all the different integer combinations of v and x, the third column has the corresponding value of W(v,x).

E.g. For v=[-1 0 1] and x=[-1 0 1], I have -

     v  x  W(v,x) 
A =  -1 -1  3 
     -1  0  4 
     -1  1  0 
      0 -1  10 
      0  0  1 
      0  1  25 
      1 -1  2 
      1  0  1 
      1  1  5

How do I obtain an approximation for P(x)=integral W(v,x)dv from this matrix in Matlab? Clearly, the integral should be replaced effectively by some summation on this matrix, but which?

W(v,x) is a phase-space density of a physical system and P(x) is the marginal distribution of the position x.

Bonus question: How do I plot W(v,x) itself with respect to both v and x (but as something like a heatmap with colours please, not using scatter3...)?


Solution

  • If I understand your question correctly, we can go this way:

    1. For getting integral: there are built-in function trapz for getting numeric integral. Now we need to create function P(x) - because we need to get trapz for current x. Here it is:

      fun = @(x) trapz( A(find(A(:,2)==x),3 ))
      

      The find function choose indexes of A for current x and then uses trapz to get integral. That is anonymous function returns integral value of W(v,x)dv for current value of x:

      fun(-1)
      ans =
         12.5000
      
    2. Second question: about plotting. We can use surf function to get 3D surface. But this function uses matrices of coordinates, not vectors, so lets reshape view of A matrix:

      v = reshape( A(:,1), [3,3])
      x = reshape( A(:,2), [3,3])
      W = reshape( A(:,3), [3,3])
      surf(v,x,W)
      

      We use [3, 3] here for this: I want to plot v as X axis, x as Y axis and W as Z. v has 3 specific values and x has 3 too. So we need to create 3x3 matrix of values! enter image description here That's not very beautiful in due to little amount of points, but I hope it is what you want.