Serving the Quantitative Finance Community

 
User avatar
Amin
Topic Author
Posts: 2416
Joined: July 14th, 2002, 3:00 am

Re: Breakthrough in the theory of stochastic differential equations and their simulation

September 25th, 2023, 10:42 pm

Inhuman, cruel and animal doggies of American army who retard thousands of innocent people without any humanity, remorse, regret, or sense of empathy because corrupt doggies want money at all cost. In a society that celebrates money (It is a quite successful model but these faults are obvious), incompetent and non-skilled people like army doggies who like many others have even a stronger desire to make money but when money cannot be made by corrupt doggies by doing anything right or pursuing any skill, they resort to finding evil ways like obliging the rich jewish billionaire and his rich cohort to make money. American army is the most corrupt segment of American society.
If some Army is professional, patriotic and civilized, people in this profession have to be praised even more than doctors, engineers and scientists. But when an army is full of inhuman, cruel and animal doggies, it is hard to appreciate anything about it. I would never have used any bad word for American army if they had any regard about respecting sanctity and dignity of so many human beings they retard without any remorse.
I demand corrupt American army and defense to stop cruel mind control of human beings all over the world and I would be the first one to praise and appreciate them. 
I have absolutely no personal animosity with anyone in American army but I want them to end their animal and cruel practices in otherwise civilized society of twenty first century.
You think life is a secret, Life is only love of flying, It has seen many ups and downs, But it likes travel more than the destination. Allama Iqbal
 
User avatar
Amin
Topic Author
Posts: 2416
Joined: July 14th, 2002, 3:00 am

Re: Breakthrough in the theory of stochastic differential equations and their simulation

September 26th, 2023, 5:57 am

When second Iraq war was started, every European country and United Nations were thoroughly convinced that Iraq had no chemical/nuclear weapons. But doggies of American army bitterly wanted to please the Jews and start the war. Every European nation tried to stop US from going to war but American army was thoroughly adamant in their claim that Saddam Hussein was making Chemical weapons.  Finally, Europeans refused to become part of any war coalition but America started the war anyway because American army bitterly wanted to attack Iraq. 
Later no chemical weapons were found. Europeans and UN were perfectly validated in their observations that Iraq had no chemical weapons plan. 
Doggie Generals of American Army started making frivolous claims that, "we had talked to a chemistry professor who had told us that Iraq was making chemical weapons. Since chemistry professor knew chemistry, we believed that Iraq must have chemical weapons." Truth was that doggies bitterly wanted to please the same Jewish Goldfather and his cohort of mind control and get material rewards.
Collin Powell who travelled to every European capital with the Dossier of evidence later demanded an inquiry was he was misled by doggies of US army who presented false evidence and his dossier had false information. This inquiry was never made. Doggies are powerful.
More than a million people were killed in the aftermath for war that doggies started for their lust for money.
Why have I written about Iraq war something which is seemingly out of context with my cause to end mind control of innocent human beings?
In order to convince people that doggie Generals of US Army are well versed in the art of lying and deception. They lie, contrive and deceive in the matters of greatest international importance where lives of millions of people are at stake as in the case of Iraq War. So lying about an ordinary individual like me is not difficult for doggies and is an everyday matter.
Since I have stepped on the tails of inhuman cruel doggies of US army, I am sure they would resort to their time tested tradition of lying and deception to thwart and fail me in my human cause. I am afraid that they can make false claims that I want to hurt people of US or people of US army or something more ingenious.
I have always told people on this forum several times that most American people are some of the best people on the planet. They are nice, kind, humane and very enterprising. It is no wonder that despite poor leadership and despite huge corruption by American army, America continues to make progress in so many dimensions and this is because Americans in general are very hard working and enterprising people and they are also well versed in all modern science and technology and all of this is possible because a large number of individuals try to excel in different spheres of life and in different professions. There is a tradition of excellence in United States.
I am also truly indebted to a large number of good Americans in the universities who tried to help me in the past and my systematic persecution decreased due to their protests. It is a matter of pride for me that some of my research work could possibly be helpful to American researchers and scientists and I would love to continue to try to help in any way I can. I know that my insight about practical applications of my research is very limited but it is really heartening to know that many people would apply them with their large and diverse understanding of practical problems in various fields. I admire all good people of United States.
Even though I openly say that America is probably the best country in the world in so many ways, I love my country since my childhood and I want to live the rest of my life in my own country. I do not want to go and live anywhere else in the world. 
But again I want friends to know that I never believe in negative patriotism. I would like to do something interesting in my country but at the same time I would love to help other people in my profession from different countries in any way I can. We are all human beings and we are all one family.
You think life is a secret, Life is only love of flying, It has seen many ups and downs, But it likes travel more than the destination. Allama Iqbal
 
User avatar
Amin
Topic Author
Posts: 2416
Joined: July 14th, 2002, 3:00 am

Re: Breakthrough in the theory of stochastic differential equations and their simulation

September 26th, 2023, 8:02 am

When I talk about American Army, I mean officers at high rank and decision makers. It is the foot soldiers and lower ranks that actually fight in the battle field and they are the ones who do no corruption at all. 
Also my observation from Pakistan Army is that larger the rank of the general, more corruption they do. And I am very sure that same applies to generals in the American army. 
By all means give great respect to army soldiers, veterans and  others who have fought for America in the actual bettlefield. They all deserve your great respect.
It is the corrupt officers and higher ranks that are actually corrupt doggies and this is the actual scum in American army. Scum never goes to actual battle field and actually swindles all the money. Being the scum in American army is one of the most lucrative professions in United States. Of course scum cannot become billionaire because of corruption but corrupt scum handily makes more money than most doctors, scientists and successful people in other leading professions. 
You think life is a secret, Life is only love of flying, It has seen many ups and downs, But it likes travel more than the destination. Allama Iqbal
 
User avatar
Amin
Topic Author
Posts: 2416
Joined: July 14th, 2002, 3:00 am

Re: Breakthrough in the theory of stochastic differential equations and their simulation

September 26th, 2023, 8:36 pm

Friends, I am back on track with my matlab programs. I am being able to run my programs for various SDE parameters. As I ran my programs again, I notice there are some robustness issues on the boundaries of the grid and I am trying to fix them. I hope there are no unforeseen problems now and I hope to post my programs again very soon possibly in a day. If there are no more special problems, I plan to continue my work as before and move on to a non-uniform grid and then towards regressions on the model.
You think life is a secret, Life is only love of flying, It has seen many ups and downs, But it likes travel more than the destination. Allama Iqbal
 
User avatar
Amin
Topic Author
Posts: 2416
Joined: July 14th, 2002, 3:00 am

Re: Breakthrough in the theory of stochastic differential equations and their simulation

September 27th, 2023, 5:51 am

Friends, my claim that mind control agents were changing my programs was wrong. I had not used matlab surf command carefully before and was mistakenly thinking that (Xn,Vm,Z(n,m)) would be the proper order. I earlier worked with both V and X both around one on a square grid so got under the very wrong impression that my program was working fine. They do sometimes alter my programs though very rarely(admittedly not about my public research) and it was in back of my mind otherwise I would have worked till the end to find the error. It was very silly.
However I am very happy that I became upset and it gave me courage to say things I have been unable to say before. I stand behind every single word I said about bad jews and their network and also what I said about American army. I own every single word I said.
I am sorry to friends that continuity of our research was disturbed but I am still relieved that I was able to speak the truth about wrong things in our society and I hope that it would help my cause of saving innocent people from cruel mind control of American army pigs. Again I want to very clearly say that everything I said about bad people is perfectly true to my knowledge and I stand behind it and I am very happy that I was able to say it.
If anything unforeseen does not happen, I am very determined to take my research towards perfect conclusion and want to work hard towards it and just like always I will continue to distribute all my matlab programs.
You think life is a secret, Life is only love of flying, It has seen many ups and downs, But it likes travel more than the destination. Allama Iqbal
 
User avatar
Amin
Topic Author
Posts: 2416
Joined: July 14th, 2002, 3:00 am

Re: Breakthrough in the theory of stochastic differential equations and their simulation

September 27th, 2023, 7:50 am

These are very initial programs for review by friends and many changes are required. I want to improve the programs in several ways. I also have some ideas about an alternative algorithm for construction of analytic density and associated 2D Hermite series. But just sharing the initial programs. 
I have given parameters for maximum values for asset and volatility and please use them freely to not let the grid become large and keep it manageable. This can be changed once a non-uniform grid is introduced. If the paths are not too large, keep the grid coarse for proper construction of numerical density.
When volatility goes smoothly into zero and a lot of mass is at zero, there is a mismatch between analytic and numerical density at zero. I am currently using CEV power .95 for volatility. I am sure this can be improved later. 

I hope that this is just start of the work and there is a lot of good research work to follow.

Here is the program.
function [] = SVMonteCarloDensityWith2DHermitesWilmott02()

%Copyright Ahsan Amin. Infiniti derivatives Technologies.
%Please fell free to connect on linkedin: linkedin.com/in/ahsan-amin-0a53334 

%Volatility and asset equations are
%Vt(0)=V0;
%dVt=(mu1*Vt^beta1 + mu2*Vt^beta2)*dt+sigma0 * Vt^gamma * dz1(t)
%Xt(0)=X0
%dXt=(a*Xt^alpha1+b*Xt^alpha2)*dt+sigma1* Vt^gammaV * Xt^gammaX *dz2(t);
%<dz1(t),dz2(t)>=rho *dt

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

X0=1.0;
alpha1=0;
alpha2=1;
kappaX=0;
thetaX=1.0;
a=kappaX*thetaX;
b=-kappaX;

rho=0;
sigma1=1;
gammaV=.5;
gammaX=.85;
%%%%%%%%`%%%%%
V0=.25;
thetaV=.06;
kappaV=1.50;
mu1=kappaV*thetaV;
mu2=-kappaV;
beta1=0;
beta2=1;
gamma=.95;
sigma0=.80;
%%%%%%%%%%%%

dt=.03125/2;
Tt=64*4;

T=Tt*dt;


%%%%%%%%%%%%%%%%%%%%%%%%5

%rng(29079137, 'twister')
rng(15898837, 'twister')
paths=500000;
V(1:paths)=V0;  %Original process monte carlo.
X(1:paths)=X0;
Random1(1:paths)=0;
Random2(1:paths)=0;

for tt=1:Tt
Random1=randn(size(Random1));
Random2=randn(size(Random2));


X(1:paths)=X(1:paths)+ ...
    (a* X(1:paths).^alpha1 + b* X(1:paths).^alpha2)* dt + ...
    sqrt(1-rho^2)* sigma1* V(1:paths).^gammaV.* X(1:paths).^gammaX .*Random1(1:paths) * sqrt(dt) + ...
    rho* sigma1* V(1:paths).^gammaV .*X(1:paths).^gammaX .*Random2(1:paths)*sqrt(dt) + ...
    (a*alpha1* X(1:paths).^(alpha1-1)+b*alpha2* X(1:paths).^(alpha2-1)).* ...
    (((a* X(1:paths).^alpha1 + b* X(1:paths).^alpha2)* dt^2/2)+ ...
    (sqrt(1-rho^2)* sigma1* V(1:paths).^gammaV.* X(1:paths).^gammaX .*Random1(1:paths) *(1-1/sqrt(3)).*dt^1.5+ ...
    rho* sigma1* V(1:paths).^gammaV .*X(1:paths).^gammaX .*Random2(1:paths)*(1-1/sqrt(3)).*dt^1.5))+ ...
    .5*(a*alpha1*(alpha1-1)* X(1:paths).^(alpha1-2)+b*alpha2*(alpha2-1).* X(1:paths).^(alpha2-2)).* ...
    ( sigma1^2* V(1:paths).^(2*gammaV).* X(1:paths).^(2*gammaX)) *dt^2/2 + ...
    sqrt(1-rho^2)* sigma1* V(1:paths).^gammaV.*gammaX.* X(1:paths).^(gammaX-1).* ...
    ((a* X(1:paths).^alpha1 + b* X(1:paths).^alpha2).*Random1(1:paths) * 1/sqrt(3).* dt^1.5 + ...
    sqrt(1-rho^2)* sigma1* V(1:paths).^gammaV.* X(1:paths).^gammaX .*(Random1(1:paths).^2-1) * dt/2 + ...
    rho* sigma1* V(1:paths).^gammaV .*X(1:paths).^gammaX .*Random1(1:paths).*Random2(1:paths)*dt/2)+ ...
    .5*sqrt(1-rho^2)* sigma1* V(1:paths).^gammaV.*gammaX.*(gammaX-1).* X(1:paths).^(gammaX-2).* ...
    (sigma1^2* V(1:paths).^(2*gammaV).* X(1:paths).^(2*gammaX) .*Random1(1:paths).*1/sqrt(3).*dt^1.5)+ ...
    sqrt(1-rho^2)* sigma1*gammaV.* V(1:paths).^(gammaV-1).* X(1:paths).^(gammaX).* ...
    ((mu1.*V(1:paths).^beta1 + mu2.*V(1:paths).^beta2).*Random1(1:paths) * 1/sqrt(3).* dt^1.5 + ...
    sigma0*V(1:paths).^gamma.*Random1(1:paths).*Random2(1:paths)*dt/2)+ ...
    .5*sqrt(1-rho^2)* sigma1*gammaV.*(gammaV-1).* V(1:paths).^(gammaV-2).* X(1:paths).^(gammaX).* ...
    (sigma0^2*V(1:paths).^(2*gamma).*Random1(1:paths)*1/sqrt(3)*dt^1.5)+ ...
    sqrt(1-rho^2)* sigma1*gammaV.* V(1:paths).^(gammaV-1).*gammaX.* X(1:paths).^(gammaX-1).* ...
    rho.* sigma1.* V(1:paths).^gammaV .*X(1:paths).^gammaX .*sigma0.*V(1:paths).^gamma.*Random1(1:paths)*1/sqrt(3)*dt^1.5+ ...
    rho* sigma1* V(1:paths).^gammaV.*gammaX.* X(1:paths).^(gammaX-1).* ...
    ((a* X(1:paths).^alpha1 + b* X(1:paths).^alpha2).*Random2(1:paths) * 1/sqrt(3).* dt^1.5 + ...
    sqrt(1-rho^2)* sigma1* V(1:paths).^gammaV.* X(1:paths).^gammaX .*Random1(1:paths).*Random2(1:paths) * dt/2 + ...
    rho* sigma1* V(1:paths).^gammaV .*X(1:paths).^gammaX .*(Random2(1:paths).^2-1)*dt/2)+ ...
    .5*rho* sigma1* V(1:paths).^gammaV.*gammaX.*(gammaX-1).* X(1:paths).^(gammaX-2).* ...
    (sigma1^2* V(1:paths).^(2*gammaV).* X(1:paths).^(2*gammaX) .*Random2(1:paths).*1/sqrt(3).*dt^1.5)+ ...
    rho* sigma1*gammaV.* V(1:paths).^(gammaV-1).* X(1:paths).^(gammaX).* ...
    ((mu1.*V(1:paths).^beta1 + mu2.*V(1:paths).^beta2).*Random2(1:paths) * 1/sqrt(3).* dt^1.5 + ...
    sigma0*V(1:paths).^gamma.*(Random2(1:paths).^2-1)*dt/2)+ ...
    .5*rho* sigma1*gammaV.*(gammaV-1).* V(1:paths).^(gammaV-2).* X(1:paths).^(gammaX).* ...
    sigma0^2.*V(1:paths).^(2*gamma).*Random2(1:paths) * 1/sqrt(3).* dt^1.5+ ...
    rho* sigma1*gammaV.* V(1:paths).^(gammaV-1).*gammaX.* X(1:paths).^(gammaX-1).* ...
    rho.* sigma1.* V(1:paths).^gammaV .*X(1:paths).^gammaX .*sigma0.*V(1:paths).^gamma.*Random2(1:paths)*1/sqrt(3)*dt^1.5;
    

    V(1:paths)=V(1:paths)+ ...
        (mu1.*V(1:paths).^beta1 + mu2.*V(1:paths).^beta2)*dt + ...
        sigma0*V(1:paths).^gamma .*Random2(1:paths)*sqrt(dt) + ...
        (mu1.*beta1*V(1:paths).^(beta1-1) + mu2.*beta2.*V(1:paths).^(beta2-1)).* ...
        ((mu1.*V(1:paths).^beta1 + mu2.*V(1:paths).^beta2)*dt^2/2 + ...
        sigma0*V(1:paths).^gamma .*Random2(1:paths)*(1-1/sqrt(3))*dt^1.5) + ...
        .5*(mu1.*beta1.*(beta1-1).*V(1:paths).^(beta1-2) + mu2.*beta2.*(beta2-1).*V(1:paths).^(beta2-2)).* ...    
        sigma0^2.*V(1:paths).^(2*gamma).*dt^2/2 + ...
        sigma0*gamma*V(1:paths).^(gamma-1) .* ...
        ((mu1.*V(1:paths).^beta1 + mu2.*V(1:paths).^beta2).*Random2(1:paths).*1/sqrt(3)*dt^1.5 + ...
        sigma0.*V(1:paths).^gamma .*(Random2(1:paths).^2-1)*dt/2) + ...
        .5*sigma0*gamma*(gamma-1).*V(1:paths).^(gamma-2) .* ...
        sigma0^2.*V(1:paths).^(2*gamma) .*Random2(1:paths).*1/sqrt(3)*dt^1.5;


end


X(X<0)=0;
V(V<0)=0;
X=real(X);
V=real(V);

SVolMeanAnalytic=thetaV+exp(-kappaV*T)*(V0-thetaV)
SVolMean=sum(V(1:paths))/paths
AssetMeanAnalytic=thetaX+exp(-kappaX*T)*(X0-thetaX)
AssetMean=sum(X(1:paths))/paths

str=input("Look at monte carlo Mean of X and V")

% NoOfBins=300;
% MaxCutOff=300;
% [XDensity0,IndexOutX0,IndexMaxX0] = MakeDensityFromSimulation_Infiniti_NEW(X,paths,NoOfBins,MaxCutOff );
% plot(IndexOutX0(1:IndexMaxX0),XDensity0(1:IndexMaxX0),'g');
% 
% str=input("Look at 1D graph of X density");

NoOfBinsV=75;
NoOfBinsX=150;
MaxCutOffX=7;
MaxCutOffV=0.65;

%[XYProbDensity2D,XYProbMass2D,XProbMass,VProbMass,IndexMaxX,IndexMaxV,Xn,Vm,Zx,Zv]= Find2DDensityAndCDFFromSimulation_Infiniti(X,V,paths,NoOfBinsX,NoOfBinsV);
[XYProbDensity2D,XYProbMass2D,XProbMass,VProbMass,IndexMaxX,IndexMaxV,Xn,Vm,Zx,Zv,BinSizeX,BinSizeV]= Find2DDensityAndCDFFromSimulation_Infiniti04(X,V,paths,NoOfBinsX,NoOfBinsV,MaxCutOffX,MaxCutOffV);
%[XYProbDensity2D,XYProbMass2D,XProbMass,YProbMass,IndexMaxX,IndexMaxY,Xn,Ym,Zx,Zy,BinSizeX,BinSizeY]
%Please look at matlab file of above function for etailed explanation of
%all its input and output parameters.



surf(Xn,Vm,transpose(XYProbDensity2D))

xlabel('Asset'), ylabel('Stochastic Volatility'), zlabel('Probability Density 2D')
title('Graph of 2D Density of Asset and SV from Parzen Density')


MeanX=sum(X(1:paths))/paths
MeanV=sum(V(1:paths))/paths

str=input("Look at 3D graphs of the numerical density")


%Below calculate Coeffs of 2D Hermite products.
%Oxn is order of Zx hermites. ranges from 1 to 8 in indexing or (0 to 7 in reality)
%Ovm is order of Zv hermites. ranges from 1 to 8 in indexing or (0 to 7 in reality)
%HermiteH(Oxn,Zx(nn,mm)) is a small function that caclulates Oxn order hermite at Zx(nn,mm) value
%HermiteH(Ovm,Zv(mm)) caclulates Ovm order hermite at Zy(mm) value
%Below Eq calculates 2D Coeffs with inner product of coordinates of 2D uniform grid with 
%product of X and V hermites.


C(1:8,1:8)=0;

for Oxn=1:8
    for Ovm=1:8
        for nn=1:(NoOfBinsX+1)
            for mm=1:(NoOfBinsV+1)
                C(Oxn,Ovm)=C(Oxn,Ovm)+1/factorial(Oxn-1)/factorial(Ovm-1).*Xn(nn).*Vm(mm).*HermiteH(Oxn-1,Zx(nn,mm)).*HermiteH(Ovm-1,Zv(mm)).*XYProbMass2D(nn,mm);
                %C(Oxn,Ovm)=C(Oxn,Ovm)+1/factorial(Oxn-1)/factorial(Ovm-1).*Xn(nn).*HermiteH(Oxn-1,Zx(nn,mm)).*HermiteH(Ovm-1,Zv(mm)).*XYProbMass2D(nn,mm);
                
            end
        end
    end
end
        
C

%If C are 2D hermite Coeffs as calculated in above equation
%Hermite Series representation of 2D random variable would be given by
%double sum equation below as 
%Sum_{n=0}^{n=7} Sum_{m=0}^{m=7} C(n+1,m+1) He(n,Zx) He(m,Zy);

str=input("Look at Coefficients of 2D Hermite products")


%Below we calculate 1D Hermite Series for volatility.
%We only need to calculate the mean (zeroth hermite) aH0 of the volatility
%random variable.


aH0=0;
aH(1:7)=0;
for nn=1:IndexMaxV
aH0  =  aH0+Vm(nn).*VProbMass(nn);    
aH(1)=aH(1)+Vm(nn).*Zv(nn).*VProbMass(nn);
aH(2)=aH(2)+.5*Vm(nn).*(Zv(nn).^2-1).*VProbMass(nn);
aH(3)=aH(3)+1/6*Vm(nn).*(Zv(nn).^3-3*Zv(nn)).*VProbMass(nn);
aH(4)=aH(4)+1/24*Vm(nn).*(Zv(nn).^4-6*Zv(nn).^2+3).*VProbMass(nn);
aH(5)=aH(5)+1/120*Vm(nn).*(Zv(nn).^5-10*Zv(nn).^3+15*Zv(nn)).*VProbMass(nn);
aH(6)=aH(6)+1/720*Vm(nn).*(Zv(nn).^6-15*Zv(nn).^4+45*Zv(nn).^2-15).*VProbMass(nn);
aH(7)=aH(7)+1/720/7*Vm(nn).*(Zv(nn).^7-21*Zv(nn).^5+105*Zv(nn).^3-105*Zv(nn)).*VProbMass(nn);
end


[a0,a] = ConvertHCoeffsToZCoeffs(aH0,aH,7);

a0
a

str=input("Look at volatility coefficients ----a0,a(1:7)");



%Below we calculate the Hermite series of asset by integrating over the
%volatility to get marginal series of asset. When we integrate over
%volatility, all volatility hermites(alone or as products with X-hermites)
%in the 2D asset series go to zero(they all have mean zero). 
%So basically all volatility terms go to zero and we are left with 
%one dimensional array containing terms from 2D series where 
%volatility hermites in the product had order zero. This order zero is
%represented by volatility index one in the 2d array of coeffs C(n,m);
%So our 1D asset series is from the line C(1:8,1) where volatilit hermites
%are zero. This is very interesting.

bH0=C(1,1);
bH(1:7)=C(2:8,1);

%Above bH0 is asset Hermite series coeff of zeroth order hermite and
%bH(1:7) represent hermites from 1st to 7th order.

bH0=bH0/aH0;
bH=bH/aH0;

%Above, we need to normalize the asset hermite series with expected value
%of volatility since we integrated over the volatility. Just like in
%Conditioning when we integrate over a variable to get marginal density and
%then we have to divide with the integral of the integrated density as
%well. Similalry, Here we have to divide with the mean/expected value of 
%volatility hermite series(which is integral of the volatility hermite series)
%to normalize our asset equation hermite coefficients.

%Below we convert from hermite series to Z-series of the asset for
%plotting.
[b0,b] = ConvertHCoeffsToZCoeffs(bH0,bH,7);

bH0
bH

b0
b

str=input("Look at 1D Coeffs of X")


[pXY2D] = Calculate2VariablesOn2DGridFrom2DZSeries(Xn,Vm,NoOfBinsX,NoOfBinsV,C,a0,a,BinSizeX,BinSizeV);

Xn0(1:NoOfBinsX-1)=Xn(2:NoOfBinsX);
Vm0(1:NoOfBinsV-1)=Vm(2:NoOfBinsV);
pXY2D0(1:NoOfBinsX-1,1:NoOfBinsV-1)=pXY2D(2:NoOfBinsX,2:NoOfBinsV);

%surf(Xn,Vm,transpose(pXY2D))
surf(Xn0,Vm0,transpose(pXY2D0))
xlabel('Asset'), ylabel('Stochastic Volatility'), zlabel('Probability Density 2D')
title('Graph of 2D Density of Asset and SV From Z-series')
str=input("Look at graph of Analytic 2D density");
XYProbDensity2D0(1:NoOfBinsX-1,1:NoOfBinsV-1)=XYProbDensity2D(2:NoOfBinsX,2:NoOfBinsV);
surf(Xn0,Vm0,transpose(pXY2D0-XYProbDensity2D0))
xlabel('Asset'), ylabel('Stochastic Volatility'), zlabel('Difference of Probability Density 2D')
title('Graph of 2D Difference of Parzen Density and Z-series Analytic Density of Asset and SV on the same grid')
str=input("Look at difference of Analytic and Numerical density graphs")




end
.
.
.
function [XYProbDensity2D,XYProbMass2D,XProbMass,YProbMass,IndexMaxX,IndexMaxY,Xn,Ym,Zx,Zy,BinSizeX,BinSizeY] = Find2DDensityAndCDFFromSimulation_Infiniti02(X,Y,Paths,NoOfBinsX,NoOfBinsY,MaxCutOffX,MaxCutOffY)
                                                                                                                                                            

%This function processes joint observations of two stochastic variables X and Y from 
%common monte carlo paths and makes a two dimensional density of X and Y.

%Explanation of Input Parameters
%X and Y are 1D arrays so that nth element of each makes a joint 
%observation of nth monte carlo path X(n) and Y(n) form joint nth observation of two stochastic variables.
%paths ---- total monte carlo paths.
%NoOfBinsX ---- Uniform Grid size in X for density construction. You can alter
%this. But you will have to make sure to alter grid sizes in main program after
%this function has returned results.
%NoOfBinsY ---- Uniform Grid size in Y for density construction. You can alter
%this.

%Explanation of output Parameters
%XYProbDensity2D --- is 2D numerically constructed Parzen density on input
%grid in X and Y.
%XYProbMass2D, ------ is integrated version of 2D parzen density that gives
%two D probability mass on input grid in X and Y.
%XProbMass  ---- is probability mass of 1D marginal density in X on NoOfBinsX Grid.
%YProbMass  ---- is probability mass of 1D marginal density in Y on NoOfBinsY Grid.
%IndexMaxX ----- This is size of grid in X and corresponds to NoOfBinsX and
%is NoOfBinsX+1. Naturally This is also size of array that returns
%X-values on the X-axis of 2D uniform grid.
%IndexMaxY ----- This is size of grid in X and corresponds to NoOfBinsY and
%is  NoOfBinsY+1. Naturally This is also size of array that returns
%Y-values on Y-axis of the 2D uniform grid.
%Xn--- is an array of size IndexMaxX that contains values of X at the
%center of each grid cell on X-axis on 2D uniform grid.
%Ym--- is an array of size IndexMaxY that contains values of Y at the
%center of each grid cell on Y-axis on 2D uniform grid.
%Zx----- this is a 2D array. This returns the value of standard normal associated with Xn for every
%Ym grid cell. Given a certain grid cell Ym, we calculate CDF along all X grid
%cells with the same Ym. This is derived from inversion of univariate CDF on X-grid
%or X-axis associated with each Y grid cell. Zx(n,m) is a 2D array.
%Zy----- this is the value of standard normal associated with Ym at center
%of each grid cell on Y-axis. This is derived from inversion of univariate CDF on Y-grid
%or Y-axis. This is 1D array. This is used for vol which is independent of
%X.


meanX=sum(X(:))/Paths;
meanY=sum(Y(:))/Paths;
sdX=sqrt(sum((X(:)-meanX).^2));
sdY=sqrt(sum((Y(:)-meanY).^2));





Xmin=99999;
Xmin1=999998;
Xmin2=999997;
Xmax=0;
Xmax1=1*10^(-12);
Xmax2=1*10^(-11);

Ymin=999999;
Ymin1=999998;
Ymin2=999997;
Ymax=0;
Ymax1=1*10^(-12);
Ymax2=1*10^(-11);
%meanX=0;
%meanY=0;
for p=1:Paths

if(isreal(X(p))==0)
    X(p)
    str=input("Look at X(p)");
end
if(isreal(Y(p))==0)
    Y(p)
    str=input("Look at Y(p)");
end    

if(X(p)>MaxCutOffX)
X(p)=MaxCutOffX;
end
if(Y(p)>MaxCutOffY)
Y(p)=MaxCutOffY;
end

if(Xmin>real(X(p)))
    Xmin2=Xmin1;
    Xmin1=Xmin;
    Xmin=real(X(p));
elseif(Xmax<real(X(p)))
    Xmax2=Xmax1;
    Xmax1=Xmax;
    Xmax=real(X(p));
end

if(Ymin>real(Y(p)))
    Ymin2=Ymin1;
    Ymin1=Ymin;
    Ymin=real(Y(p));

elseif(Ymax<real(Y(p)))
    Ymax2=Ymax1;
    Ymax1=Ymax;
    Ymax=real(Y(p));
end
end


%Xmin
%Xmin1
%Xmin2

%Xmax
%Xmax1
%Xmax2
%str=input('Look at Xmin and Xmax');
%IndexMax=NoOfBins+1;
BinSizeX=(Xmax2-Xmin2)/NoOfBinsX;
BinSizeY=(Ymax2-Ymin2)/NoOfBinsY;
%IndexMax=floor((Xmax-Xmin)/BinSize+.5)+1
IndexMaxX=floor((Xmax2-Xmin2)/BinSizeX+.5)+1;
IndexMaxY=floor((Ymax2-Ymin2)/BinSizeY+.5)+1;

BinSizeX
BinSizeY
IndexMaxX
IndexMaxY

Xmax2
Xmax1
Xmax
Xmin2
Xmin1
Xmin

Ymax2
Ymax1
Ymax
Ymin2
Ymin1
Ymin

str=input("Look at max-min parameters of numerical density")

%IndexOut(1:IndexMax)=Xmin+(0:(IndexMax-1))*BinSize;
Xn(1:IndexMaxX)=Xmin2+(0:(IndexMaxX-1))*BinSizeX;
Ym(1:IndexMaxY)=Ymin2+(0:(IndexMaxY-1))*BinSizeY;

XYProbMass2D(1:IndexMaxX,1:IndexMaxY)=0.0;
XYProbDensity2D(1:IndexMaxX,1:IndexMaxY)=0.0;
XProbMass(1:IndexMaxX)=0;
YProbMass(1:IndexMaxY)=0;

for p=1:Paths
indexX=real(floor(real(X(p)-Xmin2)/BinSizeX+.5)+1);
%indexXp=real(X(p)-Xmin2)/BinSizeX+1;


if(real(indexX)<1)
indexX=1;
end
if(real(indexX)>IndexMaxX)
indexX=IndexMaxX;
end


indexY=real(floor(real(Y(p)-Ymin2)/BinSizeY+.5)+1);
%indexYp=real(Y(p)-Ymin2)/BinSizeY;

if(real(indexY)<1)
indexY=1;
end
if(real(indexY)>IndexMaxY)
indexY=IndexMaxY;
end



%XDensity(index)=XDensity(index)+1.0/Paths/BinSize;
XProbMass(indexX)=XProbMass(indexX)+1.0/Paths;
YProbMass(indexY)=YProbMass(indexY)+1.0/Paths;
XYProbMass2D(indexX,indexY)=XYProbMass2D(indexX,indexY)+1.0/Paths;
XYProbDensity2D(indexX,indexY)=XYProbDensity2D(indexX,indexY)+1.0/Paths/BinSizeX/BinSizeY;


end

%BinSizeX
%BinSizeY

%str=input("Look at BinSizes");

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%5
%Below calculate Zx by inversion of 1D X-CDF associated with each volatility grid cell
%so that resulting array is 2D.
XProbSum(1:IndexMaxX,1:IndexMaxY)=0;
XCDF(1:IndexMaxX,1:IndexMaxY)=0;
Zx(1:IndexMaxX,1:IndexMaxY)=0;

for mm=1:IndexMaxY
    XProbSum(1,mm)=XYProbMass2D(1,mm)*.5;
    if(YProbMass(mm)>0)
        XCDF(1,mm)=XProbSum(1,mm)./YProbMass(mm);
        
    else
        XCDF(1,mm)=0;
    end
    if((XCDF(1,mm)>0) && (XCDF(1,mm)<1))
        Zx(1,mm)=norminv(XCDF(1,mm));
    elseif (XCDF(1,mm)<.000000001)
        Zx(1,mm)=-5;
    elseif (XCDF(1,mm)>.999999999)
        Zx(1,mm)=5;
    end
    for nn=2:IndexMaxX
        XProbSum(nn,mm)=XProbSum(nn-1,mm)+XYProbMass2D(nn-1,mm)*.5+XYProbMass2D(nn,mm)*.5;
        if(YProbMass(mm)>0)
            XCDF(nn,mm)=XProbSum(nn,mm)./YProbMass(mm);
        else
            XCDF(nn,mm)=0;
        end
        if((XCDF(nn,mm)>0) && (XCDF(nn,mm)<1))
            Zx(nn,mm)=norminv(XCDF(nn,mm));
        elseif (XCDF(nn,mm)<.000000001)
            Zx(nn,mm)=-5;
        elseif (XCDF(nn,mm)>.999999999)
            Zx(nn,mm)=5;
        end
        if (Zx(nn,mm)<-5)
            Zx(nn,mm)=-5;
        end
        if (Zx(nn,mm)>5)
            Zx(nn,mm)=5;
        end
    end
end
    
    

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%



YCDF(1:IndexMaxY)=0;
Zy(1:IndexMaxY)=0;
YCDF(1)=YProbMass(1)*.5;
Zy(1)=norminv(YCDF(1));
for nn=2:IndexMaxY
    YCDF(nn)=YCDF(nn-1)+YProbMass(nn-1)*.5+YProbMass(nn)*.5;
    Zy(nn)=norminv(YCDF(nn));
    if (Zy(nn)<-5)
            Zy(nn)=-5;
    end
    if (Zy(nn)>5)
            Zy(nn)=5;
    end
end

end
.
.
.
function [pXY2D] = Calculate2VariablesOn2DGridFrom2DZSeries(XX,YY,NoOfBinsX,NoOfBinsY,C,a0,a,BinSizeX,BinSizeY)


%YY are Y-coordinates of 2D Grid. In our case this is stochastic
%volatility and is univariate and independent of asset so it can be
%inverted as a univariate Z-series. its Z-series is given by a0 and a.

[Zy] = CalculateZgivenXAndZSeriesC7(YY,a0,a);


for mm=1:NoOfBinsY+1
    if(Zy(mm)>5.0)
        Zy(mm)=5.0;
    end
    if(Zy(mm)<-5.0)
        Zy(mm)=-5.0;
    end
end
Zy


str=input("Look at Zy");

bH(1:NoOfBinsY+1,1:7)=0;
bH0(1:NoOfBinsY+1)=0;
for mm=1:NoOfBinsY+1
    %bH0(mm)=C(1,1);
    for Ovm=1:8
        for Oxn=2:8
        
            
            bH(mm,Oxn-1)=bH(mm,Oxn-1)+C(Oxn,Ovm).*HermiteH(Ovm-1,Zy(mm));
            
        end
        
        bH0(mm)=bH0(mm)+C(1,Ovm).*HermiteH(Ovm-1,Zy(mm));
    end
end


for mm=1:NoOfBinsY+1
   bH0(mm)=bH0(mm)./YY(mm);
   for Oxn=1:7
   bH(mm,Oxn)=bH(mm,Oxn)./YY(mm);
   end
end

for mm=1:NoOfBinsY+1
    [b0,b]=ConvertHCoeffsToZCoeffs(bH0(mm),bH(mm,1:7),7);
    [Zx0] = CalculateZgivenXAndZSeriesC7(XX,b0,b);
    Zx(1:NoOfBinsX+1,mm)=Zx0(1:NoOfBinsX+1);
end


for mm=1:NoOfBinsY+1
    for nn=1:NoOfBinsX+1
    if(Zx(nn,mm)>5.0)
        Zx(nn,mm)=5.0;   %possible error block
    end
    if(Zx(nn,mm)<-5.0)
        Zx(nn,mm)=-5.0;   %possible error block
    end 
%    if(Zx(nn+1,mm)<Zx(nn,mm))
%        Zx(nn+1,mm)=Zx(nn,mm);   %possible error block
%    end
    end
end




% Hermite2DSeries(1:NoOfBinsX+1,1:NoOfBinsY+1)=0;
% % 
%  for nn=1:NoOfBinsX+1
%      for mm=1:NoOfBinsY+1
%      for Oxn=1:8
%      for Ovm=1:8
%            Hermite2DSeries(nn,mm)=Hermite2DSeries(nn,mm)+C(Oxn,Ovm).*HermiteH(Oxn-1,Zx(nn,mm)).*HermiteH(Ovm-1,Zy(mm));
%      end
%      end
%      end
%  end
 
 

%Below we make calculations of density in original coordinates by doing 
%change of probability derivative with respect to gaussian density.
DinvfXX(1:NoOfBinsX+1,1:NoOfBinsY+1)=0;
for mm=1:NoOfBinsY+1
    for nn=2:NoOfBinsX
    
        DinvfXX(nn,mm) = (Zx(nn + 1,mm) - Zx(nn - 1,mm))/(XX(nn + 1) - XX(nn - 1));
    
    end
    DinvfXX(NoOfBinsX+1,mm)=DinvfXX(NoOfBinsX,mm);
    DinvfXX(1,mm)=DinvfXX(2,mm);
end


DinvfYY(1:NoOfBinsY+1)=0;
for mm=2:NoOfBinsY
    
        DinvfYY(mm) = (Zy(mm+1) - Zy(mm-1))/(YY(mm + 1) - YY(mm - 1));
    
end
DinvfYY(NoOfBinsY+1)=DinvfYY(NoOfBinsY);
DinvfYY(1)=DinvfYY(2);


pXY2D(1:NoOfBinsX+1,1:NoOfBinsY+1)=0;
for mm=2:NoOfBinsY
    for nn=2:NoOfBinsX
    
    pXY2D(nn,mm) = (normpdf(Zx(nn,mm),0, 1))*abs(DinvfXX(nn,mm)).*(normpdf(Zy(mm),0, 1))*abs(DinvfYY(mm));
 
    end
end


end

.
.
.
function [Z] = CalculateZgivenXAndZSeriesC7(X,c0,c)

Z(1:length(X))=0;

for nn=1:length(X)
    %X(nn)
    %c0
        Z(nn)=CalculateZgivenXAndZSeriesBisection2C7(X(nn),c0,c);
    %Z(nn)    
end

end

.
.
.
function [Z]=CalculateZgivenXAndZSeriesBisection2C7(Y,c0,c)

%This program takes a value of Y as input and returns the corresponding 
%value of Z given the parameters of the Z-series.
%It first see if the value of Y is greater than median or smaller than
%median. If the value of Y is greater than median Z has to be greater than
%zero and if the value of Y is smaller than median Z has to be negative.
% Once we know whether Z is greater than Z or smaller than Zero, we divide
% the program into two repeated if loops, one for Z greater than zero and 
%one for Z smaller than zero. 
%If Z is positive, we have to assign it a valid band of upper and
%lower values between which it can be found. 
%This is done sequentially by calculating Y(.5) at Z=.5, if Y is smaller 
%than Y(.5), we know that true value of Y will lie between 0 and .5
%so we assign a lower band value Za=0 and upper band value Zb=.5
%If however Y is greater than Y(.5), we know that true Z should be larger 
%than .5 and we again evaluate Y(1) at Z=1, If Y is smaller than Y(1), 
%we know that true value of Z should be between .5 and 1, and we assign a 
%lower band value Za=.5 and upper band value of Zb=1. Through repeated 
%elseif we continue until we find an appropriate band value of Z.

%Once we have determined the upper band value of Zb and lower band value of
%Za, We find Zmid=(Za+Zb)/2. We then evaluate Y(Zmid) and if Y < Y(Zmid) we
%move the upper band to Zb=Zmid while keeping Za at older value. However if
% Y> Y(Zmid), we assign the lower band value equal to Zmid as Za=Zmid. In a
% few iterations, this bisection procedure gets extremely close to true
% value of Z for which Y(Z)=Y

if(Y>c0)
    %Z=.5;
    if(Y< EvaluateZSeriesC7(c0,c,.5))
        Za=0;
        Zb=.5;
    elseif(Y< EvaluateZSeriesC7(c0,c,1))
        Za=.5;
        Zb=1;
    elseif(Y< EvaluateZSeriesC7(c0,c,1.5))
        Za=1;
        Zb=1.5;  
    elseif(Y< EvaluateZSeriesC7(c0,c,2))
        Za=1.5;
        Zb=2;
    elseif(Y< EvaluateZSeriesC7(c0,c,2.5))
        Za=2;
        Zb=2.5; 
    elseif(Y< EvaluateZSeriesC7(c0,c,3))
        Za=2.5;
        Zb=3;
    elseif(Y< EvaluateZSeriesC7(c0,c,3.5))
        Za=3;
        Zb=3.5;
    elseif(Y< EvaluateZSeriesC7(c0,c,4))
        Za=3.5;
        Zb=4;
    elseif(Y< EvaluateZSeriesC7(c0,c,4.5))
        Za=4;
        Zb=4.5;
    elseif(Y< EvaluateZSeriesC7(c0,c,5))
        Za=4.5;
        Zb=5;
    elseif(Y< EvaluateZSeriesC7(c0,c,5.5))
        Za=5;
        Zb=5.5;  
    elseif(Y< EvaluateZSeriesC7(c0,c,6))
        Za=5.5;
        Zb=6;      
    else
        Za=6;
        Zb=12;
    end
    for mm=1:35
        Zmid=.5*Za+.5*Zb;
        if(Y<  EvaluateZSeriesC7(c0,c,Zmid))
            Zb=Zmid;
        elseif(Y>  EvaluateZSeriesC7(c0,c,Zmid))
            Za=Zmid;
        end
    Z=.5*Za+.5*Zb; %This can be taken out of the loop and does not 
    %necessarily have to be calculated again and again. Has to be
    %calculated at least once after the loop to determine the true value of
    %Z.
    
    end
    
end
        
if(Y<c0)
    if(Y>  EvaluateZSeriesC7(c0,c,-.5))
        Za=-.5;
        Zb=0;
    elseif(Y>  EvaluateZSeriesC7(c0,c,-1))
        Za=-1;
        Zb=-.5;
    elseif(Y>  EvaluateZSeriesC7(c0,c,-1.5))
        Za=-1.5;
        Zb=-1;  
    elseif(Y>  EvaluateZSeriesC7(c0,c,-2))
        Za=-2;
        Zb=-1.5;
    elseif(Y>  EvaluateZSeriesC7(c0,c,-2.5))
        Za=-2.5;
        Zb=-2.0; 
    elseif(Y>  EvaluateZSeriesC7(c0,c,-3))
        Za=-3;
        Zb=-2.5;
    elseif(Y>  EvaluateZSeriesC7(c0,c,-3.5))
        Za=-3.5;
        Zb=-3.0;
    elseif(Y>  EvaluateZSeriesC7(c0,c,-4))
        Za=-4;
        Zb=-3.5;
    elseif(Y>  EvaluateZSeriesC7(c0,c,-4.5))
        Za=-4.5;
        Zb=-4.0;
    elseif(Y>  EvaluateZSeriesC7(c0,c,-5))
        Za=-4.5;
        Zb=-5.0;  
    elseif(Y>  EvaluateZSeriesC7(c0,c,-5.5))
        Za=-5.0;
        Zb=-5.5; 
    elseif(Y>  EvaluateZSeriesC7(c0,c,-6.0))
        Za=-5.5;
        Zb=-6.0;     
    else
        Za=-12.0;
        Zb=-6.0;
    end
    for mm=1:35
        Zmid=.5*Za+.5*Zb;
        if(Y< EvaluateZSeriesC7(c0,c,Zmid))
            Zb=Zmid;
        elseif(Y> EvaluateZSeriesC7(c0,c,Zmid))
            Za=Zmid;
        end
    Z=.5*Za+.5*Zb; %This can be taken out of the loop and does not 
    %necessarily have to be calculated again and again. Has to be
    %calculated at least once after the loop to determine the true value of
    %Z
    end
    
end        
        
end

.
.
.
function [X] = EvaluateZSeriesC7(c0,c,Z)

X=c0+c(1)*Z+c(2)*Z.^2+c(3)*Z.^3+c(4)*Z.^4+c(5)*Z.^5+c(6)*Z.^6+c(7)*Z.^7;

end

.
.
.
This is the last few outputs you should see on your screen when you run this program.

a0 =
   0.051856848795143
a =
  Columns 1 through 3
   0.025308447061333   0.006232951335150   0.001147500450485
  Columns 4 through 6
   0.000877636380746   0.000394818257614  -0.000020803817895
  Column 7
  -0.000017810888272
Look at volatility coefficients ----a0,a(1:7)
bH0 =
   0.999230567144020
bH =
  Columns 1 through 3
   0.608033135481096   0.159292825114088   0.027749297555118
  Columns 4 through 6
   0.002484662411384  -0.001193699130006  -0.000894232267067
  Column 7
  -0.000332249594428
b0 =
   0.860805213270085
b =
  Columns 1 through 3
   0.541765963280562   0.104144398627779   0.004800081440268
  Columns 4 through 6
   0.015898146417386   0.005783542352976  -0.000894232267067
  Column 7
  -0.000332249594428
Look at 1D Coeffs of X
str =
     []
Zy =
  Columns 1 through 3
  -2.634004451734654  -2.167319894339016  -1.506064279914426
  Columns 4 through 6
  -0.874060740963614  -0.403909689433931  -0.028130067155871
  Columns 7 through 9
   0.288271345467365   0.559033179619291   0.791351162428327
  Columns 10 through 12
   0.991148471650376   1.164052401240042   1.315111681913550
  Columns 13 through 15
   1.448541791520256   1.567709385250055   1.675248979525350
  Columns 16 through 18
   1.773211661733512   1.863199918814644   1.946475835640740
  Columns 19 through 21
   2.024043532168434   2.096710658392112   2.165133997455996
  Columns 22 through 24
   2.229853354067018   2.291316897368233   2.349900274937681
  Columns 25 through 27
   2.405921167512133   2.459650482451252   2.511321048172249
  Columns 28 through 30
   2.561134433439292   2.609266346429649   2.655870948154188
  Columns 31 through 33
   2.701084328444267   2.745027330245648   2.787807863038324
  Columns 34 through 36
   2.829522812062351   2.870259626048210   2.910097647218208
  Columns 37 through 39
   2.949109233573836   2.987360712962982   3.024913200300944
  Columns 40 through 42
   3.061823303280107   3.098143736897327   3.133923863621021
  Columns 43 through 45
   3.169210173226020   3.204046713879507   3.238475484708033
  Columns 46 through 48
   3.272536798591318   3.306269623215485   3.339711907850869
  Columns 49 through 51
   3.372900903174013   3.405873481293384   3.438666464011476
  Columns 52 through 54
   3.471316967836174   3.503862775753078   3.536342747749586
  Columns 55 through 57
   3.568797284948232   3.601268866150349   3.633802681411908
  Columns 58 through 60
   3.666447395451542   3.699256085390516   3.732287414874008
  Columns 61 through 63
   3.765607132278092   3.799290020535409   3.833422487899952
  Columns 64 through 66
   3.868106089015782   3.903462431546359   3.939640211632650
  Columns 67 through 69
   3.976825640718744   4.015258517880284   4.055258214437345
  Columns 70 through 72
   4.097268263758451   4.141938997338002   4.190297468179779
  Columns 73 through 75
   4.244153386120161   4.307333394208399   4.391834105372254
  Column 76
   5.000000000000000
Look at Zy
Look at graph of Analytic 2D density
Look at difference of Analytic and Numerical density graphs

.
.
.
These are some output graphs in no particular order

Image

Image

Image

Image

Image

You think life is a secret, Life is only love of flying, It has seen many ups and downs, But it likes travel more than the destination. Allama Iqbal
 
User avatar
Amin
Topic Author
Posts: 2416
Joined: July 14th, 2002, 3:00 am

Re: Breakthrough in the theory of stochastic differential equations and their simulation

September 27th, 2023, 8:24 am

In the case of SV modelling, one SDE of volatility is independent of Asset and is fully univariate while second variable asset is bivariate and depends on both itself and volatility. Since volatility is independent, it really helps in inversion of 2D Hermite Series to find the density of asset but if both variables would be dependent on each other and both variables are bivariate, the inversion of Hermite Series would be more difficult. But I still believe that it can be handled by a 2D Bisection like algorithm. I would like to work on it but it would be slightly off topic for current project in mathematical finance. However if some friends have any such completely bivariate model, and they want to solve it, please let me know and give me the data and I would try to write the 2D bisection algorithm in my free time. I am sure it can be done.
You think life is a secret, Life is only love of flying, It has seen many ups and downs, But it likes travel more than the destination. Allama Iqbal
 
User avatar
Amin
Topic Author
Posts: 2416
Joined: July 14th, 2002, 3:00 am

Re: Breakthrough in the theory of stochastic differential equations and their simulation

September 27th, 2023, 8:56 pm

Friends, I want to very briefly explain the algorithm for inversion of 2D Hermite series to get the bivariate probability density on the specified grid of asset [$]X_n[$] and volatility [$]V_m[$] on a uniform or non-uniform rectangular grid. My notation might be different from matlab program.
We have an independent univariate Z-series/Hermite Series of volatility. 

[$]V(Z_v)\, = \, c_0 \, + c_1 \, Z_v \,+ c_2 \, {Z_v}^2 \,+ c_3 \, {Z_v}^3 \,+ c_4 \, {Z_v}^4 \,+ c_5 \, {Z_v}^5 \,+ c_6 \, {Z_v}^6 \,+ c_7 \, {Z_v}^7 \,[$]
Equation(1)

Each point on stochastic volatility grid is represented as [$]V_m[$]  and this variable has a numerical value [$]v_m[$] at that point.

We first invert the above Z-series of volatility in Eq(1) to find the value of [$]Z_v(m)[$] on each volatility grid point [$]V_m[$]. I have used the more stable bisection root finding procedure for that. So that [$]V_m(Z_v(m))\, =\, v_m[$]

Then we have equation for bivariate Hermite Series that is equal to bivariate product of stochastic volatility and asset on each cell of the rectangular grid.
This equation is given as

[$]U(Z_x,Z_y) \, =\, \sum\limits_{k_1=0}^7 \, \sum\limits_{k_2=0}^7 \, C(k_1,k_2) \, H_{k_1}(Z_x) \, H_{k_2}(Z_y) [$]  Equation(2)
We can also write the above bivariate product U on the grid points  (X_n,V_m) as
[$]U(Z_x(X_n,V_m),Z_y(V_m)) \,=X_n\, V_m [$]

We already know the value of [$]Z_v(V_m)[$] on all of the volatility grid points as we earlier did that. We substitute the numerical value of [$]Z_v(V_m)[$] in the two dimensional hermite series (sequentially for each grid point [$]V_m[$]) to get a Z_series/Hermite Series only in [$]Z_x[$]. This new series in only in [$]Z_x[$] since numerical value of [$]Z_v(V_m)[$] has been substituted in the bivariate Hermite Series. This resulting Hermite Series in [$]Z_x[$] is the product of univariate hermite Series of [$]X_n[$] , representing the X_axis associated with the particular value of stochastic volatility associated with [$]Z_v[$] that we had substituted in the 2D Hermite Series.
To get a pure Z-series of X-axis assiociated with the particular stochastic volatility point we divide the above resulting Z-series with the value of [$]V_m=v_m[$] on that volatility grid point. And after the division with volatility at that grid point, it gives us a true Z-series of X-axis associated with the particular grid point [$]V_m[$]

Writing the above explanation in the form of Equations, we have if the numerical value of stochastic volatility at mth grid point [$]V_m=v_m[$], and associated numerical value of underlying [$]Z_v=z_m[$], we substitute it in the 2D Hermite Equation U and divide by [$]v_m[$]. This division is done since U here is product of [$]X_n[$] and [$]V_m[$] and we want to divide by [$]V_m=v_m[$] to get the equation of only [$]X_n[$] out of the equation of the product.

[$]X(Z_x(n,m)) \, =\Large \, \frac{\, U(Z_x(X_n,V_m),Z_v(V_m)) \, \LARGE {|}^{Z_v(V_m)=z_m}}{v_m}[$]  Equation(3)

The above formula gives us a univariate Hermite Series for that line of asset grid that is associated with mth volatility grid point [$]V_m[$]. This univariate Z-series can easily be inverted with our bisection algorithm to find [$]Z_x[$] associated with [$]X_n[$] line on mth volatility grid point.

The above formula is repeated for each of the points [$]V_m[$] on the stochastic volatility grid and it gives us a Z-series that is valid on the X points associated only for mth volatility grid cell. 
if total number of volatility grid points is M, We can then invert each of the M dependent Z-series to get the value of [$]Z_x(n,m)[$]. Each of the Z-series is valid for X points associated with that particular volatility grid point that was substituted in the Equation(3) to get the value of Z-series. All of these M Z-series are inverted with the bisection root-finding algorithm that I have previously posted by using a loop that has M iterations.

So with above analytics, we calculate [$]Z_x(n,m)[$] and [$]Z_v(m)[$] associated with each grid point [$]X_n[$] and [$]V_m[$] on the rectangular grid.

After this we calculate the derivatives [$]\frac{dZ_v(m)}{dV_m}[$] and [$]\frac{dZ_x(n,m)}{dX_n}[$] on all of the grid points using simple finit differences. These derivatives are then used to calculate the Jacobian of the density of bivariate variable [$]X_nV_m[$] with respect to two standard normals. 
Since stochastic volatility did not have any dependence on [$]Z_x[$], our jacobian is simply product of above two derivatives.

We get the Equation of Bivariate Density as

[$]p(X_n,V_m)\, = \, p_N(Z_x(n,m)) \, \frac{dZ_x(n,m)}{dX_n} \,  \, p_N(Z_y(m)) \, \frac{dZ_y(m)}{dV_m} \,[$] 

Here in the RHS of Equation [$]p_N[$] denotes standard normal density.
You think life is a secret, Life is only love of flying, It has seen many ups and downs, But it likes travel more than the destination. Allama Iqbal
 
User avatar
Amin
Topic Author
Posts: 2416
Joined: July 14th, 2002, 3:00 am

Re: Breakthrough in the theory of stochastic differential equations and their simulation

September 28th, 2023, 8:28 am

Friends, I have another idea to find coefficients of 2D Hermite Series for SV models by least squares polynomial regression. I had first tried this idea on univariate data with limited success. But later when I standardized the univariate data and then ran the regression, my results were very good. I was just using several tens of data points. In monte carlo setting, I believe that we can only use several tens of thousands of paths to get a very stable estimate of 2D Hermite series coefficients with least squares polynomial regression method. In my last program, I was using half a million paths for reasonably stable numerical density but I really think that in a polynomial regression setting, we could do well only with fifty thousand or even less paths or even smaller if we use low discrepancy sequences. 
I really want to try this idea within next two days and post the matlab program and results of my work on the forum.
You think life is a secret, Life is only love of flying, It has seen many ups and downs, But it likes travel more than the destination. Allama Iqbal
 
User avatar
Amin
Topic Author
Posts: 2416
Joined: July 14th, 2002, 3:00 am

Re: Breakthrough in the theory of stochastic differential equations and their simulation

September 29th, 2023, 7:54 am

Friends, I really believe that in a civilized society of twenty first century there is absolutely no place for connivance and machinations of a group of rich jews led by jewish billionaire Godfather of mind control who have retarded, at whim, thousands of talented intelligent people in America and abroad that they did not like. These cruel people are inhuman animals even though they have learnt to adjust extremely well within American society. 
These corrupt crooks are very powerful and they have already complete control over American army and a large number of politicians. Due to this huge influence and power, it is very difficult to do anything to undo the cruel mind control network led by jewish billionaire and his rich jewish cohort.
I have suggested several times previously to good American people if nothing else works, please just ask some patriotic American investigative journalists in major newspapers to run a comprehensive story about mind control network of jewish billionaire and the nexus between rich jews and American defense and about tens of thousands of people that have been retarded over past thirty years.
I want to invite patriotic American journalists from CNN, Washington Post, New York Times, Guardian or some other large, reputed and responsible US or European media outlets to run a thorough story on machinations of evil jews led by jewish billionaire and the  nexus of bad jews with American defense and their cruel mind control network that has retarded tens of thousands of very intelligent Americans and foreigners all over the world.
Since good people may not have the power to forcefully end the evil practices of the powerful nexus between evil Jewish billionaire and US army crooks, only way to end the cruel and subhuman practices of the nexus seems to be exposing them openly in US public and the public of other countries all across the world. 
I think it is important for us to make sure there is a better future of American children and children all across the world and we leave them a legacy in which fairness and justice prevails as opposed to machinations, trickery and deceit that has been unleashed in the society by the jewish billionaire and his rich jewish cohort.
You think life is a secret, Life is only love of flying, It has seen many ups and downs, But it likes travel more than the destination. Allama Iqbal
 
User avatar
Amin
Topic Author
Posts: 2416
Joined: July 14th, 2002, 3:00 am

Re: Breakthrough in the theory of stochastic differential equations and their simulation

September 29th, 2023, 11:06 am

Friends, In my last matlab program with finding 2D Density of SV model with hermite polynomials, I used Zx(n,m) i.e there is a different Zx for each volatility grid cell. Though it seemed to work decently, I think I need to calculate a single universal Zx for all of the volatility grid cells. I want to try this and see how the new program works with a single universal Zx. 
For calculation of 2D hermite series of SV model with polynomial regression that I suggested yesterday, I have to take a single universal Zx. I have worked out details of that but I want to see how a single Zx goes with the previous program.
You think life is a secret, Life is only love of flying, It has seen many ups and downs, But it likes travel more than the destination. Allama Iqbal
 
User avatar
Amin
Topic Author
Posts: 2416
Joined: July 14th, 2002, 3:00 am

Re: Breakthrough in the theory of stochastic differential equations and their simulation

September 29th, 2023, 6:22 pm

Mind control technology was developed with the help of large number of neurotransmitter experts and some of them at the top of their profession. And a large number of neurotransmitter experts are provided data of neurotransmitter interactions in the brain of several mind control victims (especially when victims are hard to retard) and they advise American defense, after analyzing the data, which neurotransmitters to block and how to reconnect the victims's brain in order to successfully retard the victim. These neurotransmitter experts continue to analyze the data to make suggestions to mind control agencies. These neurotransmitter experts are extremely exorbitantly paid by American defense for their services. 
And sometimes when a more well known victim has to be retarded and people around the target are absolutely opposed to let any body harm the target victim, American defense tells the responsible people around the target victim that a certain leading neurotransmitter expert has told them with perfect certainty that the target victim has evil neurotransmitters. And since very exorbitant amounts of money are paid to neurotransmitter experts, these neurotransmitter experts make a solemn face and blurt totally complete lies in the most professional way to tell the people around the victim that the victim indeed has extremely evil neurotransmitters. These dishonest neurotransmitter experts tell complete lies of the sort that whoever they knew with certain neurotransmitter always turned evil. Since non-professional people have no way to know, many of them can be deceived easily by professional dishonesty of the experts. 
Has anything like this ever happened to someone around you? You have been thoroughly wronged by dishonest experts.
You think life is a secret, Life is only love of flying, It has seen many ups and downs, But it likes travel more than the destination. Allama Iqbal
 
User avatar
Amin
Topic Author
Posts: 2416
Joined: July 14th, 2002, 3:00 am

Re: Breakthrough in the theory of stochastic differential equations and their simulation

September 29th, 2023, 8:48 pm

Friends, I really think that we need to do more work on representation of two random variables that are dependent or those that are correlated with appropriate Z-series.
When one of the random variables is dependent on other, this dependence can be captured by making coefficients of its Z-series a function of the other random variable. For example if independent Z-series of volatility is given as

[$]\, V \, = \, b_0 \, + \, b_1 \, Z_v \,+ \, b_2 \, {Z_v}^2 \,+ \, b_3 \, {Z_v}^3 \,+ \, b_4 \, {Z_v}^4 \,+ \, b_5 \, {Z_v}^5 \,+ \, b_6 \, {Z_v}^6 \,+ \, b_7 \, {Z_v}^7 \,[$] 

Then dependent Z-Series of asset can be represented as

[$]\, X \, =\, c_0 \, + \, c_1(Z_v) \, Z_x \,+ \, c_2(Z_v) \, {Z_x}^2 \,+ \, c_3(Z_v) \, {Z_x}^3 \,+ \, c_4(Z_v) \, {Z_x}^4 \,+ \, c_5(Z_v) \, {Z_x}^5 \,+ \, c_6(Z_v) \, {Z_x}^6 \,+ \, c_7(Z_v) \, {Z_x}^7 \,[$]

So when a Z-series is dependent on other, this dependence has to be captured by making its coefficients a function of the other Z-series variable. These coefficients have to be found so that dependence structure is properly captured.

However when a random variable is non-dependent but correlated with another random variable, this correlation has to be captured by 

[$]\, X \, =\, c_0 \, + \, c_1 \, Z_x \,+ \, c_2\, {Z_x}^2 \,+ \, c_3 \, {Z_x}^3 \,+ \, c_4 \, {Z_x}^4 \,+ \, c_5 \, {Z_x}^5 \,+ \, c_6 \, {Z_x}^6 \,+ \, c_7 \, {Z_x}^7 \, \\
+\, \, d_1 \, Z_y \,+ \, d_2\, {Z_y}^2 \,+ \, d_3 \, {Z_y}^3 \,+ \, d_4 \, {Z_y}^4 \,+ \, d_5 \, {Z_y}^5 \,+ \, d_6 \, {Z_y}^6 \,+ \, d_7 \, {Z_y}^7[$]

Coefficients in above Z-series are chosen so that marginal moments and correlations are simultaneously retrieved.

So dependence goes into coefficients and is captured by making coefficients of dependent Z-series function of the other Z-series, while correlation is captured by appending one Z-series with other Z-series whose coefficients are chosen to capture the correlation structure.

When a random variable is both dependent and correlated, I think we could possibly represent it as a Z-series as

[$]\, X \, =\, c_0 \, + \, c_1(Z_y) \, Z_x \,+ \, c_2(Z_y)\, {Z_x}^2 \,+ \, c_3(Z_y) \, {Z_x}^3 \,+ \, c_4(Z_y) \, {Z_x}^4 \,+ \, c_5(Z_y) \, {Z_x}^5 \,+ \, c_6(Z_y) \, {Z_x}^6 \,+ \, c_7(Z_y) \, {Z_x}^7 \, \\
+\, \, d_1(Z_y) \, Z_y \,+ \, d_2(Z_y)\, {Z_y}^2 \,+ \, d_3(Z_y) \, {Z_y}^3 \,+ \, d_4(Z_y) \, {Z_y}^4 \,+ \, d_5(Z_y) \, {Z_y}^5 \,+ \, d_6(Z_y) \, {Z_y}^6 \,+ \, d_7(Z_y) \, {Z_y}^7[$]

This would be very interesting how to empirically calculate the coefficients of Z-series of dependent random variable with a certain dependence structure and same for correlated Z-series and dependent correlated Z-series. 
I will try to think more about it but if we can somehow find numerical algorithms to construct dependent and correlated Z-series, it will be a great boost for our research towards practical applications of Z-series. I really think this would be really remarkable if we can construct algorithms for calculation of dependent and correlated Z-series.

There are other very interesting things we can do with further research. For example in order to represent a Z-series that encompasses a reasonably large class of SDEs, we can possibly make coefficients of the Z-series a function of the SDE parameters so that when parameters of the SDE change, coefficients of the Z-series are automatically altered(as they are functions of parameters) to capture the density of the SDE with particular parameters so that a single general Z-series could be used to represent an SDE variable with SDE parameters over a reasonably large range. This can be probably be impossible for a large time step/frame but I think it can be done for an arbitrary initial distribution (with a given Z-series) and an arbitrary SDE(with a reasonably large range of parameters) over a small time step of 1/16 or so (for most of the SDEs we use in finance) so that coefficients of resulting Z-series are functions of the SDE parameters.

And especially if we can do that with Z-series of stochastic volatility and dependent asset in an SV model.
You think life is a secret, Life is only love of flying, It has seen many ups and downs, But it likes travel more than the destination. Allama Iqbal
 
User avatar
Amin
Topic Author
Posts: 2416
Joined: July 14th, 2002, 3:00 am

Re: Breakthrough in the theory of stochastic differential equations and their simulation

September 30th, 2023, 10:17 am

Friends, when I first worked on hermite regressions, I really thought that they could replace AI in several important problems related to random processes (Well, AI can always absorb Mathematics and start doing more interesting things).  
Last night I continued to think of alternative and better ways to find a Z-series representation of both asset and volatility in coupled SV models given the numerical density. And then I realize that this is one of the places where we can easily use hermite regressions to find the bivariate coefficients of asset both in correlated and non-correlated models. 
I want to briefly share my ideas with friends about how to use the hermite regressions in our problem.
First we construct a bivariate grid and find the univariate Z-series of independent stochastic volatility. Then we construct a conditional Z-series of asset for every volatility grid point. There will be as many asset Z-series as points on the volatility axis of the bivariate grid, each for every volatility point.
All of the asset Z-series will have different coefficients and even median of every asset Z-series would be different. 
We will calculate moments of each of the coefficients of asset Z-series across all volatility Z-series and use our Z-series construction from moments method to fins the Z-series representation of each of the coefficients of asset Z-series. For eight coefficients, this will give us eight Z-series representing distribution of each coefficient across all of the asset Z-series(that are in turn dependent on volatility value of each volatility grid point).
Then we  hermite regress Z-series of each coefficient on Z-series of the independent volatility to extract the dependence structure of asset Z-series coefficients on independent volatility. This regression will be repeated for each of the eight asset Z-series coefficients. 
Briefly, for hermite regressions, we will calculate value of hermite polynomials of each coefficient on all of the volatility grids and pair them with values of hermite coefficients on independent volatility on all of the associated volatility grids and then run hermite regression using all these pairs.
For example this method will give us a bivariate Z-series of asset and we will automatically get conditional Z-series of asset at any volatility grid point simply by substituting the value of volatility Zv associated with that volatility grid point.
This method would work well for both correlated and non-correlated SV SDEs.
I assure friends that it would work very well and will be a wonderful application of hermite regressions method we discovered earlier this year. And I hope we will start using hermite regressions in several similar problems more freely.

Though slightly prudent due to pressure from good people, mind control agents remain totally malicious and continue to try to reconnect my brain with a sense of urgency. Many times they use another different gas in my room and other places. When inhaled this gas cause uneasiness in the stomach and also slight pain and uneasiness in the spine. The weather is becoming more pleasant and I sit in the garage/porch of our house at night since it is more open, airy  and colder but still has walls and a gate that can let heavy gas accumulate there. They continue to release gas in the garage. Even if I do not turn on the ceiling fan to avoid the gas, they still release gas that continues to accumulate in the garage. For past two nights, I completely open the gates of the house/garage to not let the gas accumulate there and spread out. When there is more gas and I am trying to think of some idea, I simply walk on the local street outside my house to be able to think better.

I have an injection due on monday. I could not get any good injection from older manufacturing dates that I knew had good injections so this time I just bought an injection with most recent manufacturing date I could get. I do not know how it would go. But I have some good coffee and some good tea that would be very helpful. I do have a feeling that mind control agents are bitterly waiting for the new injection to reconnect my brain, try to retard me somewhat and not let me work properly. But I am also determined to continue my research.
You think life is a secret, Life is only love of flying, It has seen many ups and downs, But it likes travel more than the destination. Allama Iqbal
 
User avatar
Amin
Topic Author
Posts: 2416
Joined: July 14th, 2002, 3:00 am

Re: Breakthrough in the theory of stochastic differential equations and their simulation

October 1st, 2023, 8:09 pm

Friends, I could not login to wilmott from my computer. When I try to login from main wilmott site using my computer and it says that my credentials are invalid. I could login from forums login page but it automatically logs me out when I hit any button.
I could however login from my cell phone and write this post.
You think life is a secret, Life is only love of flying, It has seen many ups and downs, But it likes travel more than the destination. Allama Iqbal