CODEWITHSUNDEEP
matlabsymbolic-mathcalculusprobability-density
Vincent.W.
Vincent.W.

Reputation: 113

Dirac probability measure in Matlab

I want to implement densities of probability measures in Matlab. For that I define density as a function handle such that the integral of some function f (given as a function handle) on the interval [a,b] can be computed by 

syms x
int(f(x)*density(x),x,a,b)

When it comes to the Dirac measure the problem is that 

int(dirac(x),x,0,b)

delivers the value 1/2 instead of 1 for all b>0. However if I type 

int(dirac(x),x,a,b)

where a<0 and b>0 the returned value is 1 as it should be. For this reason multiplying by 2 will not suffice as I want my density to be valid for all intervals [a,b]. I also dont want to distinct cases before integrating so that the code remains valid for a large class of densities.

Does someone know how I can implement the Dirac probability measure (as defined here) in Matlab?

Upvotes: 0

Views: 359

Answers (1)

Brick
Brick

Reputation: 4262

There's no unique, accepted definition for int(delta, 0, b). The problem here isn't that you're getting the "wrong" answer as that you want to impose a different convention somehow on what the delta function is than what was provided by Matlab. (Their choice is defensible but not unique.) If you evaluate this in Wolfram Alpha, for example, it will give you theta(0) - which is not defined as anything in particular. Here theta is the Heaviside function. If you want to impose your own convention, implement your own delta function.

EDIT

I see you wrote a comment on the question, while I was writing this answer, so.... Keep in mind that the Dirac measure or Dirac delta function is not a function at all. The problem you're having, along with what's described below, all relate to trying to give a functional form to something that is inherently not a function. What you are doing is not well-defined in the framework that you have in Matlab.

END OF EDIT

To put the point about conventions in context, the delta function can be defined by different properties. One is that int(delta(x) f(x), a, b) = f(0) when a < 0 < b. This doesn't tell you anything about the integral that you want. Another, which probably leads to an answer as you're getting from Matlab, is to define it as a limit. One (but not the only choice) is the limit of the zero-mean Gaussian as the variance goes to 0.

If you want to use a convention int(delta(x) f(x), a, b) = f(0) when a <= 0 < b, that probably won't get you in much trouble, but just keep in mind that it's a convention you've chosen more than a "right" or "wrong" answer relative to what you got from Matlab.

As related note, there is a similar choice to be made on the step function (Heaviside function) at x=0. There are conventions in which is is (a) undefined, (b) -1, (c) +1, and (d) 1/2. None is "wrong." This probably corresponds roughly to the choice on the Dirac function since the Heaviside is (roughly) the integral of the Dirac.

Upvotes: 1

Related Questions