• C++ Programming for Financial Engineering
    Highly recommended by thousands of MFE students. Covers essential C++ topics with applications to financial engineering. Learn more Join!
    Python for Finance with Intro to Data Science
    Gain practical understanding of Python to read, understand, and write professional Python code for your first day on the job. Learn more Join!
    An Intuition-Based Options Primer for FE
    Ideal for entry level positions interviews and graduate studies, specializing in options trading arbitrage and options valuation models. Learn more Join!

Problem using FFT to get prob dist in Matlab

Joined
5/1/09
Messages
1
Points
11
I am trying to use the Fast Fourier Transform to derive the probability distribution for an affine jump diffusion model. For some reason, I am having an enormously difficult time doing this in Matlab. I have posted my code below for trying to derive the Gaussian distribution from the Gaussian characteristic function (sort of like the Moment Generating Function). However I need to divide the resultant pdf by the interval size and still I get too much probability in the extremes (I keep getting probability's of 0.025 out for x's > 5 and x's < -5). Can someone please help me figure out what I am doing wrong.....Thanks.

In addition I have attached an excel file that contains the FFT resultant pdf and the correct pdf as generated by Matlab's 'pdf' function.

l = -10;
u = 10;
T = 512;
dx = (u-l)/T;
n = [0:T-1];
x = (l)+n*dx; %my pdf x's
t = [0:1:T-1];
s_t = 2*pi*t/(T*dx); %my s's for the Gaussian characteristic function
mu = 0;
sigma = 1;
lambda = 2;

norm = exp(s_t.*i.*mu - 0.5.*(sigma.*s_t).^2 - i.*s_t.*l); %Gaussian characteristic function (MGF except with 'it' substituted for 't'
pdfnorm = real(fft(norm))./(u-l);
plot(x,pdfnorm); %the whole distribution is wrong....it has too much probability at the extremes...and it does this even when I increase T or increase l and u...HELP!!!


This is what the pdf and psi function should look like

xlr=10; np=512; dx=2*xlr/np;
x=dx*[0:np-1]-xlr; % OR, equivalently
x=linspace(-xlr,xlr-dx,np);

normp1 = pdf('norm',x,0,1);
fnormp1 = fft(normp1);
normp1after = ifft(fnormp1);
plot (x,normp1)
 
Its been a while since Ive done this things but I believe that your problem is that you missing fftshift after the fft.
After you do fft the zero-frequency component isnt in the middle of the spectrum and you need to move it.
Im not sure if thats your problem, but from my signal processing background I can tell that this looks fishy.
After fft you have to do fftshift and after ifft you have to do ifftshift.
Try that.
 
Back
Top