Software in Mathematica for forecasting stock market, polls and insurance claims, practically anything of a stochastic nature, statistically with confidence intervals.
See greenparty.ca/blog/169 for further info.
Initialization
(* All Rights Reserved - CopyRight 2012, Richard John Belshaw, AGU,AAAS,ASA,IEEE,ACS,ACM*)
LaunchKernels[]
ClearAll[matlabFind1]
matlabFind1[lst_List,Op_:{Unequal,GreaterEqual,LessEqual},elt_:0]:Module[{x},Flatten[Position[lst,x_;Op[x,elt]]]];
ClearAll[Find2] (* new definition from 1/9/2012 )
Find2[lst_List,i_]:=Module[{L,j,count},
L=Length[lst]; count=0; Do[ If[lst[[j]]i,count=count+1],{j,1,L}];count];
Needs["Calendar`"]
getData
(
All Rights Reserved - CopyRight 2012, Richard John Belshaw, AGU,AAAS,ASA,IEEE,ACS,ACM*)
ClearAll[getData];
getData[data_List,numpts_Integer:8192,extendflag_Integer,debugflag_Integer]; ((extendflag
0||extendflag
1)&&(debugflag
1||debugflag
0)):=Module[{xa, xb, MaxRatio,ratios, difflist,phaselist, extenlst, intlist, data1,pdflist, outlst, sdifflist},
MaxRatio = Max[data]Min[data];
ratios = Ratios[data];
If[debugflag
1,Print["3a"]];
{xa, xb, extenlst, intlist, pdflist, outlst}=getPDF5[ratios,numpts,MaxRatio,extendflag,debugflag];
If[debugflag
1,Print["3b"]];
{xa, xb, ratios, extenlst, intlist, pdflist, outlst}]
getPDF5
(* All Rights Reserved - CopyRight 2012, Richard John Belshaw, AGU,AAAS,ASA,IEEE,ACS,ACM*)
ClearAll[getPDF5]; (* new definition from 3/5/2013 )
Compile-> {
getPDF5[temp_List,N2_Integer,MaxRatio_,extendflag_Integer:0,debugflag_Integer:0]
;(N2>Length[temp]&&(extendflag
0||extendflag
1)&&(debugflag
1||debugflag
0)):=Module[{sortlst,N1,xa,xb,extenlst,extenlst1,i,h0,intlist,pdflist,outlst,h7,coeffs,FlagRatioa2,FlagRatiob1,j,x0,sumval,m,mat,y,x,z,a,b,x1,x2,c,d,a5,b5,c5,d5,a7,b7,c7,d7,nlm1,a2,f7,c2,d2,f1,nlm2,b1,f8,c1,d1,z1,f2,count,temp1,n,delta,f11,f12,k,flag1,down,up,n1,n2,N3,intpdflist,pdflistfn,a3,b2,temp2,func1,pdf,int,f,temp3,z2,k1,pdflist2,h01,nlm1pars,nlm2pars},
If [debugflag
1, Print["Start GetPDF5"]];
sortlst=Sort[temp];N1=Length[temp];
xa=sortlst[[1]];xb=sortlst[[-1]];extenlst1=sortlst;extenlst
{};
temp1={};
i=1;
While[i<= N1, While[i<= N1 && (extenlst1[[i]]=0||extenlst1[[i]]
Infinity||!NumberQ[extenlst1[[i]]]),i+=1];If[i<
N1,AppendTo[extenlst,extenlst1[[i]]]];i+1];
i=2;
extenlst = Sort[extenlst];
temp1 = {};
While[i<
N1, While [i<N1 &&(extenlst[[i]]
extenlst[[i-1]]),i+=1];If[i<
N1+1,AppendTo[temp1,extenlst[[i-1]]]];i=i+1];
If [ i== N1+1 && extenlst[[i-2]]!= extenlst[[i-1]], AppendTo[temp1,extenlst[[i-1]]]];
temp1=Sort[temp1];
N1=Length[extenlst];xa=extenlst[[1]];xb=extenlst[[-1]];
(*Print[xa," , ",xb];)
If[extendflag =1,
a=xa;
b=xb;
If[debugflag
1,Print["4a"];
Print[temp1[[1;;7]]]];
nlm1pars[a2,f7,c2,d2,x] = Log[Abs[a2]] + c2*Log[x+2.+f7]+d2*Log[x+2.+f7]^2;
Off[NonlinearModelFit::nrlnum,NonlinearModelFit::cvmit,NonlinearModelFit::sszero];
nlm1 = NonlinearModelFit[Log[temp1[[1;;7]]],nlm1pars[a2,f7,c2,d2,x],{a2,f7,c2,d2},x,MaxIterations-> 30000];
z1 = nlm1["BestFitParameters"];
If[debugflag
1,Print[z1]];
a2 =Abs[z1[[1,2]]];
If[debugflag
1,Print["a2
",a2]];
FlagRatioa2=0;
If[a2<1/(50*MaxRatio),a2=1/(2*MaxRatio);FlagRatioa2=1;If[debugflag=1,Print["a2new
",a2]]];
f7 z1[[2,2]];
c2 = z1[[3,2]];
d2 = z1[[4,2]];
If [debugflag
1,Print["x1
",x1]];
If[x1<1E-20, x1=1E-20];
If[debugflag=1,Print["4b"];
Print[temp1[[-7;;-1]]]];
nlm2pars[b1,f8,c1,d1,x] = Log[Abs[b1]]+ c1*Log[x+2.+f8]-d1*Log[x+2.+f8]^2;
nlm2
NonlinearModelFit[Log[Reverse[temp1[[-7;;-1]]]],nlm2pars[b1,f8,c1,d1,x],{b1,f8,c1,d1},x,MaxIterations->30000];
(Print[z1];)
z2 = nlm2["BestFitParameters"];
b1 = Abs[z2[[1,2]]];
If[debugflag=1,Print["z2 = ",z2];Print["b1 = ",b1]];
FlagRatiob1=0;
If [b1 > 50*MaxRatio,b1 = 2*MaxRatio;FlagRatiob1=1;If[debugflag
1,Print["b1new = ",b1]]];
f8 =z2[[2,2]];
c1 = z2[[3,2]];
d1 = z2[[4,2]];
If [debugflag
1,Print["4c"]];
delta =1;
up = 0;
down = up - delta;
flag1=0;
count=0;
If [debugflag
1,Print["a2 = ",a2]];
While [((Re[Exp[ nlm1[down]]]< Re[Exp[nlm1[up]]]) && (Re[Exp[nlm1[down]]] =
Exp[nlm1[down]]) &&Re[Exp[nlm1[down ]]]< a&&Re[Exp[nlm1[down]]]>= a2&& count < 200 && FlagRatioa2=0),
If[debugflag
1,Print[Re[Exp[nlm1[down ]]]]];
up = down;
count = count +1;
down = down - delta;
];
If[debugflag
1,Print["4d "]];
up = N[up];
temp2 = {};
i
0;
While [(i>= up &&Exp[nlm1[i]]=Re[Exp[nlm1[i]]]&& Re[Exp[nlm1[i]]]=a2&& FlagRatioa20),
AppendTo[temp2,Exp[nlm1[i]]];
If[debugflag
1,Print[Exp[nlm1[i]]]];
i=i-delta;
];
If[FlagRatioa2
0,
If[temp2
{}, a3 = a, a3 = Sort[temp2][[1]]];
If [a2>0&& a3>a2,
AppendTo[temp2,a2];
If[debugflag
1,Print[a2]];
];
];
extenlst = Join[temp2,extenlst];
extenlst = N[Sort[extenlst]];
xa = extenlst[[1]],
xa = a];
If [extendflag
1,
If [debugflag
1,
Print["4e"]];
delta = 1;
up = 0;
down = up - delta;
If[debugflag
1,
Print["b1
",b1]];
count 0;
While [((Re[Exp[ nlm2[down]]]> Re[Exp[nlm2[up]]]) && (Re[Exp[nlm2[down]]] =
Exp[nlm2[down]])&&Re[Exp[nlm2[down ]]]> b&&Re[Exp[nlm2[down]]]<= b1&& count < 200 && FlagRatiob1=0),
If [debugflag
1,Print[Re[Exp[nlm2[down ]]]]];
up = down;
count = count +1;
down = down - delta;
];
If[debugflag
1,Print["4f"]];
up = N[up];
temp2 = {};
i=0;
While [(i>
up && Exp[nlm2[i]]=Re[Exp[nlm2[i]]]&&Re[Exp[nlm2[i]]]>b &&Re[Exp[nlm2[i]]]<
b1&& FlagRatiob1=0),
If [debugflag
1,Print[Exp[nlm2[i]]]];
AppendTo[temp2,Exp[nlm2[i]]];
i=i-delta;
];
If [FlagRatiob1
0,
If [ temp2 =
{}, b2 = b, b2 = Sort[temp2][[-1]]];
If [ b1 > b2,
If[debugflag=1,Print[b1]];
AppendTo[temp2,b1];
];
];
extenlst = Join[extenlst,temp2];
extenlst =N[ Sort[extenlst]];
xb = extenlst[[-1]],
xb=b];
If[debugflag
1,Print["4g"]];
N1 = Length[extenlst];
f = ConstantArray[{},N2];
func1 = ConstantArray[{},N1];
outlst = ConstantArray[0,N2];
pdflist = ConstantArray[0,N2];
xa=extenlst[[1]];
xb=extenlst[[-1]];
Clear[f2];
Do [ func1[[i]] = {i ,(SetPrecision[ extenlst[[i]],30])},{i,1,N1}];
y = Interpolation[func1];
h0 = N[(N1-1)(N2-1)];
Do [f[[i+1]] = {N[h0*i+1] ,N[Evaluate[y[h0*i+1]]]}, {i,0,N2-1}];
f2 = Interpolation[f];
temp2 = Derivative[1][f2];
Do [If[temp2[h0*i+1]
0,pdflist[[i+1]] = 0,pdflist[[i+1]]=1/N[Evaluate[temp2[h0*i+1]]]],{i,0,N2-1}];
Do [ outlst[[i+1]] = N[y[(i*h0+1)]],{i,0,N2-1}];
(* z = Fourier[pdflist];
z[[Floor[N2/2*1/5];;N2/2+1+Floor[N2/2*4/5]]]=0;
pdflist = Re[InverseFourier[z]];)
Do[If[pdflist[[i]]<0,pdflist[[i]] = 0],{i,1,N2}];
{sumval,intlist} = getIntegrate2[pdflist[[1;;N2]],1,1,N2];
pdflist = pdflist/sumval;
intlist = intlist/intlist[[-1]];
If[debugflag
1,Print["End GetPDF5"]];
{xa,xb,extenlst,intlist,pdflist,outlst}]}
Compile->{Null}
getIntegrate2
(
All Rights Reserved - CopyRight 2012, Richard John Belshaw, AGU,AAAS,ASA,IEEE,ACS,ACM*)
k = Boltzmann's constant or 1, a and b limits of integration on x axis from first to last point,eg .95 to 1.1
ClearAll[getIntegrate2];
Compile -> {
getIntegrate2[datapts_?(VectorQ[#,NumberQ]&),k_?NumberQ,a_?NumberQ,b_?NumberQ]:=Module[{N1=Length[datapts],coeffs
{14,32,12,32},h1,intlist,sumval,j,temp1},h1=((b-a)(N1-1)2/45);
intlist={7*datapts[[1]]^k*h1};sumval=intlist[[1]];
j=1;While[j<= N1-2,temp1=coeffs[[j+1-Floor[j/4]*4]]*datapts[[j+1]]^k*h1;sumval=sumval+temp1; intlist=Append[intlist,sumval];j+1];sumval=sumval+7*datapts[[N1]]^k*h1;
intlist=Append[intlist,sumval];
N[{sumval,intlist}]]}
Compile->{Null}
getIntegrate2[Range[0.05,0.95,0.05],1,1,10]
{4.52,{0.00777778,0.0788889,0.118889,0.261111,0.338889,0.552222,0.645556,0.93,1.07,1.42556,1.57222,1.99889,2.20111,2.69889,2.89889,3.46778,3.73222,4.37222,4.52}}
getDist
(
All Rights Reserved - CopyRight 2012, Richard John Belshaw, AGU,AAAS,ASA,IEEE,ACS,ACM*)
ClearAll[getDist];(new definition from 2/14/2013)getDist[lst_?(VectorQ[#,Positive[#]&]&)]:=Module[{L,i,j,k,lst1,newlist,sdifflist,sumval},lst1=lst;
L=Length[lst1];
While[j<=L,If[lst1[[j]]
0||lst1[[j]]!=NumberQ,i=j;While[i<=L,lst1[[i]]=lst1[[i+1]];i=i+1;];L=L-1,i=j;] j=j+1;];
sdifflist=Sort[Abs[Differences[Transpose[SortBy[Transpose[{Range[1,Length[lst1]],lst1}],Last]][[1]]]]];
L=L-1;
newlist=ConstantArray[0,L];
j=1;
i=1;
While[i<=L&&j<=L,k=0;
While[j<=L&&sdifflist[[j]]
i,k=k+1;
j=j+1];
newlist[[i]]=k;
i=i+1];
sumval=Total[newlist];
newlist=N[newlist/sumval]];
vec=RandomVariate[ParetoDistribution[2,2,1],15]
{6.37928,6.76598,1.25079,7.17549,1.05218,10.9294,1.69027,1.02719,1.56484,6.6485,4.03571,2.04785,1.36641,2.87586,1.00453}
getDist[vec](* value from 2/14/2012 )
{0.0714286,0.214286,0.142857,0.,0.142857,0.0714286,0.0714286,0.0714286,0.,0.0714286,0.0714286,0.0714286,0.,0.}
getPercentErrorCubic
( All Rights Reserved - CopyRight 2012, Richard John Belshaw, AGU,AAAS,ASA,IEEE,ACS,ACM*)
ClearAll[getPercentErrorCubic] (* definition from Feb 7, 2013
where TempData is original stock prices, dataarray is either Null matrix {{}} or temperatures data )
Compile -> {
getPercentErrorCubic[TempData_List,dataarray_List?MatrixQ]:=Module[{i,j,k,a=1,b,N1,N2,coeffs,vol6,vol7,L1,L,Temp,Soy,Sxoy,Sx2oy,Sx3oy,Soy2,Sxoy2,Sx2oy2,Sx3oy2,Sx4oy2,Sx5oy2,Sx6oy2,vol,volatility,vol2},{a,b}=Dimensions[dataarray];If[b
0,N1=0;N2=Length[TempData],N1=Length[TempData];N2=Length[TempData]+b];L=N2-N1; L1=N1+1;
coeffs=ConstantArray[0,4];vol6=ConstantArray[0,L];vol7=0;
Do[ If[a>1,Temp=SetPrecision[dataarray[[i,1;;L]],60],Temp=SetPrecision[TempData,60]];
Do[ If[Temp[[k]]
0,Temp[[k]] = 10^(-20)],{k,1,L}];
Soy=Total[1/Temp[[1;;L]]];
Sxoy=Total[Range[L1,N2]Temp[[1;;L]]];
Sx2oy=Total[Range[L1,N2]^2/Temp[[1;;L]]];
Sx3oy=Total[Range[L1,N2]^3/Temp[[1;;L]]];
Soy2=Total[1/Temp[[1;;L]]^2];
Sxoy2=Total[Range[L1,N2]/Temp[[1;;L]]^2];
Sx2oy2=Total[Range[L1,N2]^2/Temp[[1;;L]]^2];
Sx3oy2=Total[Range[L1,N2]^3/Temp[[1;;L]]^2];
Sx4oy2=Total[Range[L1,N2]^4/Temp[[1;;L]]^2];
Sx5oy2=Total[Range[L1,N2]^5/Temp[[1;;L]]^2];
Sx6oy2=Total[Range[L1,N2]^6/Temp[[1;;L]]^2];
coeffs[[1]]
(Sx4oy2*Sxoy2*Sx3oy2*Sxoy+Soy2*Sx3oy2*Sx4oy2*Sx2oy-Soy2*Sx3oy2^2*Sx3oy-Sx5oy2*Sx3oy2*Sxoy2*Soy+ Sx3oy2^3*Soy-Sxoy2^2*Sx4oy2*Sx3oy+2*Sx2oy2*Sxoy2*Sx3oy2*Sx3oy+Sxoy2^2*Sx5oy2*Sx2oy-Sx2oy2*Sxoy2*Sx4oy2*Sx2oy- 2*Sx2oy2*Sx4oy2*Sx3oy2*Soy-Sx2oy2*Sx5oy2*Sxoy2*Sxoy+Sx2oy2^2*Sx3oy2*Sx2oy+Soy2*Sx2oy2*Sx4oy2*Sx3oy-Soy2*Sx2oy2*Sx5oy2*Sx2oy- Sx2oy2*Sx3oy2^2*Sxoy-Sx2oy2^3*Sx3oy+Soy2*Sx5oy2*Sx3oy2*Sxoy-Soy2*Sx4oy2^2*Sxoy+Sx2oy2^2*Sx4oy2*Sxoy+Sx4oy2^2*Sxoy2*Soy- Sx3oy2^2*Sxoy2*Sx2oy+Sx2oy2^2*Sx5oy2*Soy)(2*Sx2oy2^2*Sx5oy2*Sx3oy2-Sx2oy2^3*Sx6oy2-2*Sx5oy2*Sx3oy2^2*Sxoy2-3*Sx2oy2*Sx4oy2*Sx3oy2^2+Sx2oy2^2*Sx4oy2^2+2*Sx4oy2^2*Sxoy2*Sx3oy2-Soy2*Sx2oy2*Sx5oy2^2-Soy2*Sx4oy2^3- Sxoy2^2*Sx4oy2*Sx6oy2+Sx3oy2^4-Soy2*Sx3oy2^2*Sx6oy2+Sxoy2^2*Sx5oy2^2-2*Sx2oy2*Sx5oy2*Sxoy2*Sx4oy2+2*Sx2oy2*Sxoy2*Sx3oy2*Sx6oy2+Soy2*Sx2oy2*Sx4oy2*Sx6oy2+2*Soy2*Sx3oy2*Sx4oy2*Sx5oy2);
coeffs[[2]]-(Sx3oy2*Sx5oy2*Sxoy*Sxoy2-Sx3oy2^3*Sxoy-Sx3oy*Sx5oy2*Sxoy2^2+Sx3oy*Sx3oy2^2*Sxoy2-Sx2oy2*Sx6oy2*Sxoy*Sxoy2+ Sx2oy2*Sx3oy2*Sx4oy2*Sxoy+Sx2oy2*Sx3oy*Sx4oy2*Sxoy2-Sx2oy2^2*Sx3oy*Sx3oy2+Sx2oy*Sx6oy2*Sxoy2^2-2*Sx2oy*Sx3oy2*Sx4oy2*Sxoy2+ Sx2oy*Sx2oy2*Sx3oy2^2-Soy2*Sx4oy2*Sx5oy2*Sxoy+Soy2*Sx3oy2*Sx6oy2*Sxoy-Soy2*Sx3oy*Sx3oy2*Sx4oy2+Soy2*Sx2oy2*Sx3oy*Sx5oy2+ Soy2*Sx2oy*Sx4oy2^2-Soy2*Sx2oy*Sx2oy2*Sx6oy2+Soy*Sx4oy2*Sx5oy2*Sxoy2-Soy*Sx3oy2*Sx6oy2*Sxoy2+Soy*Sx3oy2^2*Sx4oy2- Soy*Sx2oy2*Sx4oy2^2-Soy*Sx2oy2*Sx3oy2*Sx5oy2+Soy*Sx2oy2^2*Sx6oy2)(2*Sx2oy2^2*Sx5oy2*Sx3oy2-Sx2oy2^3*Sx6oy2- 2*Sx5oy2*Sx3oy2^2*Sxoy2-3*Sx2oy2*Sx4oy2*Sx3oy2^2+Sx2oy2^2*Sx4oy2^2+2*Sx4oy2^2*Sxoy2*Sx3oy2-Soy2*Sx2oy2*Sx5oy2^2-Soy2*Sx4oy2^3- Sxoy2^2*Sx4oy2*Sx6oy2+Sx3oy2^4-Soy2*Sx3oy2^2*Sx6oy2+Sxoy2^2*Sx5oy2^2-2*Sx2oy2*Sx5oy2*Sxoy2*Sx4oy2+2*Sx2oy2*Sxoy2*Sx3oy2*Sx6oy2+Soy2*Sx2oy2*Sx4oy2*Sx6oy2+2*Soy2*Sx3oy2*Sx4oy2*Sx5oy2);
coeffs[[3]]
(Soy2*Sx3oy2*Sx3oy*Sx5oy2+Soy2*Sx4oy2*Sx5oy2*Sx2oy-Soy2*Sx4oy2^2*Sx3oy-Soy2*Sxoy*Sx5oy2^2-Soy2*Sx3oy2*Sx2oy*Sx6oy2+ Soy2*Sxoy*Sx4oy2*Sx6oy2+Sx2oy2^2*Sx4oy2*Sx3oy-Sx2oy2^2*Sxoy*Sx6oy2+Sx2oy2*Sxoy2*Sx2oy*Sx6oy2-Sx2oy2*Sxoy2*Sx3oy*Sx5oy2- Sx2oy2*Sx4oy2*Sx5oy2*Soy+Sx2oy2*Sx3oy2*Soy*Sx6oy2+2*Sx2oy2*Sx5oy2*Sx3oy2*Sxoy-Sx3oy2^2*Sx2oy2*Sx3oy- Sx3oy2*Sx4oy2*Sx2oy2*Sx2oy-Sx4oy2*Sx3oy2^2*Sxoy+Sxoy2*Soy*Sx5oy2^2-Sx3oy2^2*Soy*Sx5oy2-Sxoy2*Sx3oy2*Sx5oy2*Sx2oy- Sxoy2*Soy*Sx4oy2*Sx6oy2+Sx4oy2^2*Sx3oy2*Soy+Sx3oy2^3*Sx2oy+Sx4oy2*Sxoy2*Sx3oy2*Sx3oy)(2*Sx2oy2^2*Sx5oy2*Sx3oy2- Sx2oy2^3*Sx6oy2-2*Sx5oy2*Sx3oy2^2*Sxoy2-3*Sx2oy2*Sx4oy2*Sx3oy2^2+Sx2oy2^2*Sx4oy2^2+2*Sx4oy2^2*Sxoy2*Sx3oy2- Soy2*Sx2oy2*Sx5oy2^2-Soy2*Sx4oy2^3-Sxoy2^2*Sx4oy2*Sx6oy2+Sx3oy2^4-Soy2*Sx3oy2^2*Sx6oy2+Sxoy2^2*Sx5oy2^2- 2*Sx2oy2*Sx5oy2*Sxoy2*Sx4oy2+2*Sx2oy2*Sxoy2*Sx3oy2*Sx6oy2+Soy2*Sx2oy2*Sx4oy2*Sx6oy2+2*Soy2*Sx3oy2*Sx4oy2*Sx5oy2);
coeffs[[4]](Sxoy2*Sx4oy2^2*Sx3oy+Sxoy2*Sxoy*Sx5oy2^2-Soy*Sx2oy2*Sx5oy2^2-Soy*Sx3oy2^2*Sx6oy2-Sx3oy2^2*Sx4oy2*Sx2oy- Sx5oy2*Sx3oy2^2*Sxoy+Sx3oy2*Sx4oy2^2*Sxoy+Sx2oy2^2*Sx3oy*Sx5oy2+Sx2oy2*Sx2oy*Sx4oy2^2-Sx2oy*Sx2oy2^2*Sx6oy2- 2*Sx3oy2*Sx2oy2*Sx4oy2*Sx3oy+Sx3oy2*Sx2oy2*Sx5oy2*Sx2oy-Sx2oy2*Sx4oy2*Sx5oy2*Sxoy+Sx2oy2*Sx3oy2*Sx6oy2*Sxoy- Sxoy2*Sx3oy2*Sx3oy*Sx5oy2-Sxoy2*Sx4oy2*Sx5oy2*Sx2oy+Sxoy2*Sx3oy2*Sx2oy*Sx6oy2-Sxoy2*Sxoy*Sx4oy2*Sx6oy2+ Soy*Sx2oy2*Sx4oy2*Sx6oy2+2*Soy*Sx3oy2*Sx4oy2*Sx5oy2+Sx3oy2^3*Sx3oy-Soy*Sx4oy2^3)(2*Sx2oy2^2*Sx5oy2*Sx3oy2- Sx2oy2^3*Sx6oy2-2*Sx5oy2*Sx3oy2^2*Sxoy2-3*Sx2oy2*Sx4oy2*Sx3oy2^2+Sx2oy2^2*Sx4oy2^2+2*Sx4oy2^2*Sxoy2*Sx3oy2- Soy2*Sx2oy2*Sx5oy2^2-Soy2*Sx4oy2^3-Sxoy2^2*Sx4oy2*Sx6oy2+Sx3oy2^4-Soy2*Sx3oy2^2*Sx6oy2+Sxoy2^2*Sx5oy2^2- 2*Sx2oy2*Sx5oy2*Sxoy2*Sx4oy2+2*Sx2oy2*Sxoy2*Sx3oy2*Sx6oy2+Soy2*Sx2oy2*Sx4oy2*Sx6oy2+2*Soy2*Sx3oy2*Sx4oy2*Sx5oy2);
vol=ConstantArray[0,L];
Do[vol[[j]]
(Temp[[j]]-coeffs[[1]](N1+j)^3-coeffs[[2]](N1+j)^2-coeffs[[3]](N1+j)-coeffs[[4]])^2;
vol6[[j]]vol6[[j]]+vol[[j]],{j,1,L}];
vol7=vol7+Sqrt[Total[vol[[1;;L]]]
(L-1)],{i,1,a}];
If[a>1, vol2=vol7/a;volatility=vol6/a,volatility=vol7;vol2=0];
{coeffs,volatility,vol2}]}
Compile->{Null}
generateMagPhase
(* All Rights Reserved - CopyRight 2012, Richard John Belshaw, AGU,AAAS,ASA,IEEE,ACS,ACM*)
ClearAll[generateMagPhase];(new definition from Mar 4,2013)
Compile-> {generateMagPhase[ratios_List,intlist_List,pdflist_List,outlst_List,xa_?NumberQ,xb_?NumberQ,np_Integer,steps_Integer,ppast_,debugflag_Integer];((0<=ppast<=1)&&(debugflag
1||debugflag
0)):=Module[{phasedistpdf,ft,sigma,psi,mu,x,nlm,h,h1,dB,pdB,phaselistpdf2,phaselistcdf,sumval,N3,i,j,k,x0,indxvals,g1list,g2list,rv,k1,p2,futureflag,mag,p1,phasedistcdf,indices,temp,phasedistpdf2,k2,p11,x2,N5,weights,p3},
If [debugflag
1,Print["Start GenerateMagPhase"]];
phasedistpdf=getDist[SetPrecision[ratios,20]];
If[debugflag
1,Print["*1"]];
{sumval,phasedistcdf}=getIntegrate2[phasedistpdf,1,1,Length[phasedistpdf]];
phasedistcdf=phasedistcdf/sumval;
If[debugflag
1,Print["
**2"]];
indices=matlabFind1[phasedistpdf,Unequal,0];
If[debugflag
1,Print["
**3"]];
ft[sigma_,psi_,mu_,x_]:=Re[1/sigma(1+psi*(x-mu)/sigma)^((-1/sigma)-1)];
(ft[sigma_,psi_,mu_,x_]:=Re[PDF[ParetoDistribution[sigma,psi,mu],x]];)
nlm=NonlinearModelFit[phasedistpdf,ft[sigma,psi,mu,x],{sigma,psi,mu},x,Method->"NMinimize",MaxIterations->1000,PrecisionGoal->20];
{sigma,psi,mu}=Transpose[nlm["ParameterTableEntries"]][[1]];
(* h[sigma_,psi_,mu_,x_]:=If[1+psi*(x-mu)/sigma>0,Re[1/sigma*(1+psi*(x-mu)/sigma)^((-1/sigma)-1)],0]; )
(*weights = ConstantArray[1,9*Length[phasedistpdf]];
*)
weights = ft[sigma,psi,mu,#]&
@Range[Length[phasedistpdf]+1,9*Length[phasedistpdf]];
phaselistpdf2=Join[phasedistpdf,weights];
{sumval,phaselistcdf}=getIntegrate2[phaselistpdf2,1,1,Length[phaselistpdf2]];
phaselistcdf=phaselistcdf/sumval;
N3=Length[intlist];indxvals=ConstantArray[0,Length[intlist]];
i=1;k=1;
While[i<=N3,x0
(i-1)(N3-1);
While[k<N3&&intlist[[k]]1,indxvals[[i]]=k-1,indxvals[[i]]=1.];
i=i+1];
g1list=ConstantArray[0,Length[phaselistpdf2]];
i=1;k=1;
While[i<=Length[phaselistpdf2],x0
(i-1)(Length[phaselistpdf2]-1);
While[k<Length[phaselistcdf]&&phaselistcdf[[k]]1,g1list[[i]]
(k-1),g1list[[i]]1.];
i=i+1];
g2list=ConstantArray[0,Length[indices]];
i=1;k=1;
While[i<=Length[indices],x0
(i-1)(Length[phasedistpdf]-1);
While[k<Length[indices]&&phasedistcdf[[indices[[k]]]]1,g2list[[i]]
(k-1),g2list[[i]]1.];
i=i+1];
dB=ConstantArray[0,{np,steps}];
pdB=ConstantArray[0,{np,steps}];
N5=Length[pdflist];
Do[rv=RandomReal[];
futureflag=0;
If[rv=j,mag=ratios[[Length[ratios]-(k-j)]],mag=dB[[i,j-k]]],k1=Floor[rv
(Length[phaselistpdf2]-1)],
If[k1<=Length[phasedistpdf],rv=RandomReal[];
k1=Floor[rv*(Length[indices]-1)+1.5];
k=g2list[[k1]];
p2=phasedistpdf[[indices[[k]]]],
k2=g1list[[k1]];
p2=phaselistpdf2[[k2]];
futureflag=1;rv=RandomReal[];
k1=Floor[rv*(Length[intlist]-1)+1.5];
k2=indxvals[[k1]];mag=outlst[[k2]]]];
If[magxb,p1=0,
temp
(mag-xa)(xb-xa)(N5-1)+1;
k1=Floor[temp];
x2=temp-k1;
p11=pdflist[[k1]];
If[k1{Null}
generateTemperatures
( All Rights Reserved - CopyRight 2012, Richard John Belshaw, AGU,AAAS,ASA,IEEE,ACS,ACM*)
ClearAll[generateTemperatures]; (* 11/6/2012 and 12/10/2012 and 12/11/2012 where s0 is initial price )
Compile -> {
generateTemperatures[dB_List?MatrixQ,s0_,np_Integer,steps_Integer];(Dimensions[dB]={np,steps}):=Module[{temperatures,i,j},
temperatures=ConstantArray[0,{np,steps}];
i=1;
While[i<=np,
j=2; temperatures[[i,1]]=s0*dB[[i,1]];
While[j<=steps,temperatures[[i,j]]=temperatures[[i,j-1]]*dB[[i,j]];j=j+1];
i=i+1];
temperatures]}
Compile->{Null}
generateTemperatures[RandomReal[{0.95,1.05},{10,5}],21.2,10,5]
{{20.5374,21.184,20.7382,21.2536,20.266},{21.6525,20.8237,19.9736,20.1481,19.225},{20.4032,20.3139,21.2922,20.8431,20.0757},{21.3645,20.3421,19.4224,19.8336,19.4298},{20.8459,20.998,21.9209,22.7808,23.2908},{21.2042,21.6289,21.9336,21.2845,20.7828},{22.2567,22.1804,21.8878,22.7006,23.7557},{21.1944,20.4451,21.2166,21.5168,21.7687},{20.1886,20.5268,20.2092,19.9364,19.4733},{21.7436,21.1628,21.2765,21.1087,21.144}}
calculateConfLevels
(
All Rights Reserved - CopyRight 2012, Richard John Belshaw, AGU,AAAS,ASA,IEEE,ACS,ACM*)
ClearAll[calculateConfLevels]; (* 11/11/2012 )
Compile-> {
calculateConfLevels[temperatures_List?MatrixQ,pdB_List?MatrixQ,s0_, negFlag_Integer] /; (negFlag1||negFlag
0):=Module[{np,steps,conf01,conf10,conf25,conf60,conf50,conf75,conf90,conf99,temp1,columnindx,temper2,i,j,ratio,temp2,sumval1,k,sumval2,confarray,conford,sumval,indxlist,p1,p2,columnord},{np,steps}=Dimensions[temperatures];
columnindx = ConstantArray[0,np];
conf01=ConstantArray[0,steps];
conf10=conf25=conf50=conf75=conf90=conf99=conf01;
confarray = ConstantArray[0,{np,steps}];
conford = ConstantArray[0,{np,steps}];
Do[temp1=temperatures[[1;;np,j]];columnindx=Sort[temp1];columnord = Ordering[temp1];
i=Floor[.01*np+0.5];If[i<1,i=1];
conf01[[j]]=columnindx[[i]];
i=Floor[.1*np+0.5];
conf10[[j]]=columnindx[[i]];
i=Floor[.25*np+0.5];
conf25[[j]]=columnindx[[i]];
i=Floor[.5*np+0.5];
conf50[[j]]=columnindx[[i]];
i=Floor[.75*np+0.5];
conf75[[j]]=columnindx[[i]];
i=Floor[.9*np+0.5];
conf90[[j]]=columnindx[[i]];
i=Floor[.99*np+0.5];
conf99[[j]]=columnindx[[i]];
confarray[[1;;np,j]] = columnindx;
conford[[1;;np,j]] = columnord;
,{j,1,steps}];
ratio =indxlist
ConstantArray[0,198];
p2 = ConstantArray[0,198];
p1 = ConstantArray[0,steps];
Do[p2[[j]] 1-( .5+.0025*j),{j,1,198}];
Do[k=Floor[(.5+.0025*j)*np+.5];indxlist[[j]] = p2[[j]];Do[p1[[i]]
pdB[[conford[[k,i]],i]],{i,1,steps}]; temp1 p1[[1;;steps]]*p2[[j]]
( confarray[[k,1;;steps]]-s0);
{sumval1,temp2} = getIntegrate2[temp1,1,1,steps];sumval1=sumval1/Total[p1]/p2[[j]]/steps;
k=Floor[(.5-.0025*j)np+.5];temp1 = p2[[j]]*p1[[1;;steps]](s0-confarray[[k,1;;steps]]);
Clear[temp2];
{sumval2,temp2} = getIntegrate2[temp1,1,1,steps];sumval2 = sumval2/Total[p1]/p2[[j]]/steps;
If[sumval2
0,sumval2=10^(-20)];
Clear[temp2];
temp2 = sumval1/sumval2;
ratio[[j]]=If [temp2<0,If [negFlag
1,-N[1
(1+Abs[temp2])100.]-100,N[Abs[temp2](1+Abs[temp2])*100.+100]],If[negFlag1,-N[1
(1+temp2)*100],N[temp2/(1+temp2)*100.]]],{j,1,198}];
Clear[temp2];
{sumval,temp2} = getIntegrate2[ratio,1,1,198];
sumval = sumval/198;
{conf01,conf10,conf25,conf50,conf75,conf90,conf99,sumval,indxlist,ratio,ratio[[66]],ratio[[100]],ratio[[160]],ratio[[198]]}]}
Compile->{Null}
calculateConfLevels[generateTemperatures[RandomReal[{0.95,1.05},{10,5}],21.2,10,5]]
{{20.2496,19.938,19.6515,20.2559,19.5276},{20.2496,19.938,19.6515,20.2559,19.5276},{20.8274,20.8171,20.6641,21.4057,20.8806},{21.3592,21.2926,21.535,21.5731,21.0701},{21.8998,21.7119,22.0443,21.7467,22.5259},{22.0064,21.8179,22.6244,21.8547,22.558},{22.0959,22.1243,22.889,22.9716,23.0679},{{20.2496,19.938,19.6515,20.2559,19.5276},{20.6404,20.6132,20.4042,20.476,20.5985},{20.8274,20.8171,20.6641,21.4057,20.8806},{20.9243,20.9085,21.4206,21.4383,20.884},{21.3592,21.2926,21.535,21.5731,21.0701},{21.6572,21.3579,21.7078,21.6363,21.4264},{21.6969,21.4158,22.0109,21.6639,22.3367},{21.8998,21.7119,22.0443,21.7467,22.5259},{22.0064,21.8179,22.6244,21.8547,22.558},{22.0959,22.1243,22.889,22.9716,23.0679}}}
calculateFutureValues
( All Rights Reserved - CopyRight 2012, Richard John Belshaw, AGU,AAAS,ASA,IEEE,ACS,ACM*)
ClearAll[calculateFutureValues];(s0 initial price,k7 strike price,r int rate,d discount rate,infl inflation rate,T is time ratio,{np,steps} is dimensions of Temperatures,dB,pdB, Feb 21 later 2013 *)
Compile-> {calculateFutureValues[dB_List?MatrixQ,pdB_List?MatrixQ,futurearma_List,s0_?NumberQ,k7_?NumberQ,r_?NumberQ,d_?NumberQ,infl_?NumberQ,T_?NumberQ,debugflag_Integer];((Dimensions[pdB]Dimensions[dB])&&(Length[futurearma]
Dimensions[dB][[2]])&&(debugflag
1||debugflag
0)):=Module[{np,steps,dis,confcall10,confcall50,confcall60,confcall90,confput10,confput50,confput60,confput90,Tvalforwardup,Tvalforwarddown,valueforwardup,valueforwarddown,temperatures,meanpath,dt,i,sumprob,L,h7=1.,j,newpassedthresholdforwarddown,newpassedthresholdforwardup,v,valsup,lengthsup,temp,valsdown,lengthsdown,nup,ndown,indx,count1,count2,TvalForwarddown,TvalForwardup,futurearma1,pi=4 ArcTan[1],p,phase,newpassedthresholdsindown,newpassedthresholdsinup,indxup,indxdown,valuesindownd,valuesinupd,temp1,k,callup,calldown,confcall,confput,callup1,calldown1,temp2,p1,p2,sumprob1,t1,t2,valsupord,valsdownord,sumval,callup2,calldown2},{np,steps}=Dimensions[dB];dt=T/steps;
If[debugflag
1,Print["StartCalcFutureValues"]];
temperatures=generateTemperatures[dB,s0,np,steps];
i=1;
If [debugflag
1,Print["1"]];
sumprob=ConstantArray[0,steps];L=steps;
Do[sumprob[[j]]=Re[Total[pdB[[1;;np,j]]]];If[sumprob[[j]]
= 0, sumprob[[j]]1],{j,1,steps}];
meanpath=ConstantArray[0,steps];
futurearma1=ConstantArray[0,steps];
Do[meanpath[[j]]=Total[pdB[[1;;np,j]]*temperatures[[1;;np,j]]]/sumprob[[j]];
futurearma1[[j]]=futurearma[[j]]+meanpath[[j]],{j,1,steps}];
newpassedthresholdforwarddown=ConstantArray[0,{np,steps}];
newpassedthresholdforwardup=ConstantArray[0,{np,steps}];
newpassedthresholdsindown=ConstantArray[0,steps];
newpassedthresholdsinup=ConstantArray[0,steps];
Do[v=futurearma1[[j]]-k7;
If[v>0,newpassedthresholdsinup[[j]]=v,newpassedthresholdsinup[[j]]=0];
v
-v;
If[v>0,newpassedthresholdsindown[[j]]v,newpassedthresholdsindown[[j]]=0];
Do[v=temperatures[[i,j]]-k7;
If[v>0,newpassedthresholdforwardup[[i,j]]=v,newpassedthresholdforwardup[[i,j]]=0];
v
-v;
If[v>0,newpassedthresholdforwarddown[[i,j]]v,newpassedthresholdforwarddown[[i,j]]=0],{i,1,np}],{j,1,steps}];
valueforwarddown=ConstantArray[0,steps];
valueforwardup=ConstantArray[0,steps];
valuesindownd=ConstantArray[0,steps];
valuesinupd=ConstantArray[0,steps];
j=2;dis=ConstantArray[0,steps];(*discount rate
)dis[[1]]
(1+(r+infl-d)dt);
While[j<steps,dis[[j]]=dis[[j-1]]*dis[[1]];
j=j+1];
valsup=ConstantArray[0,{np,steps}];
valsupord = ConstantArray[0,{np,steps}];
indxup=ConstantArray[0,{np,steps}];
temp1 = ConstantArray[0,np];
lengthsup=ConstantArray[0,steps];
Do[k=1;
Do[If[newpassedthresholdforwardup[[i,j]]>0,indxup[[k,j]]=i;
k=k+1],{i,1,np}];
k=k-1;
temp=Sort[newpassedthresholdforwardup[[indxup[[1;;k,j]],j]]];
valsupord[[1;;k,j]]=Ordering[newpassedthresholdforwardup[[indxup[[1;;k,j]],j]]];
lengthsup[[j]]=k;
valsup[[1;;k,j]]=temp,{j,1,steps}];
valsdown=ConstantArray[0,{np,steps}];
valsdownord = ConstantArray[0,{np,steps}];
lengthsdown=ConstantArray[0,steps];
indxdown = ConstantArray[0,{np,steps}];
Do[k=1;
Do[If[newpassedthresholdforwarddown[[i,j]]>0,indxdown[[k,j]]=i;
k=k+1],{i,1,np}];
k=k-1;
temp=Sort[newpassedthresholdforwarddown[[indxdown[[1;;k,j]],j]]];
valsdownord[[1;;k,j]]=Ordering[newpassedthresholdforwarddown[[indxdown[[1;;k,j]],j]]];
lengthsdown[[j]]=k;
valsdown[[1;;k,j]]=temp,{j,1,steps}];
confcall10=confcall50=confcall60=confcall90=confput10=confput50=confput60=confput90=ConstantArray[0,steps];
Do[nup=lengthsup[[j]];
If[nup>0,indx=Floor[0.1*nup+.5];If[indx<1,indx=1];
confcall10[[j]]=valsup[[indx,j]];
indx = Floor[.5*nup+.5];
confcall50[[j]] = valsup[[indx,j]];
indx = Floor[.6*nup+.5];
confcall60[[j]]=valsup[[indx,j]];
indx=Floor[0.9*nup+0.5];
confcall90[[j]]=valsup[[indx,j]]];
ndown=lengthsdown[[j]];
If[ndown>0,indx=Floor[0.1*ndown+.5];If [indx < 1, indx =1];
confput10[[j]]=valsdown[[indx,j]];
indx = Floor[.5*ndown+.5];
confput50[[j]]=valsdown[[indx,j]];
indx = Floor[.6*ndown+.5];
confput60[[j]] = valsdown[[indx,j]];
indx=Floor[0.9*ndown+0.5];
confput90[[j]]=valsdown[[indx,j]]],{j,1,steps}];
callup = ConstantArray[0,{198*2+1,steps}];
callup2 = ConstantArray[0,{198*2+1,steps}];
callup1 = ConstantArray[0,steps];
calldown = ConstantArray[0,{198*2+1,steps}];
calldown2 = ConstantArray[0,{198*2+1,steps}];
calldown1 = ConstantArray[0,steps];
t1 = ConstantArray[0,steps];
t2 = ConstantArray[0,steps];
p2 = ConstantArray[0,2*198+1];
Do[p
.005+.0025*k;If[p>.5,p2[[k]]1-p,p2[[k]]=p],{k,0,2*198}];
Do[nup = lengthsup[[j]];ndown=lengthsdown[[j]];Do[ p
.005+.0025*k;
If[nup>0,indx Floor[p*nup+.5];If[indx<1,indx=1];p1=pdB[[valsupord[[indx,j]],j]];callup2[[k+1,j]] = valsup[[indx,j]];callup[[k+1,j]] =p2[[k]]*p1*valsup[[indx,j]],callup[[k+1,j]]=callup2[[k+1,j]] = 0];
t1[[j]]=t1[[j]]+p1*p2[[k]];
If[ndown>0,indx=Floor[p*ndown+.5];If[indx<1,indx=1];p1=pdB[[valsdownord[[indx,j]],j]];calldown2[[k+1,j]]=valsdown[[indx,j]];calldown[[k+1,j]] = p2[[k]]*p1*valsdown[[indx,j]],calldown2[[k+1,j]]=calldown[[k+1,j]]=0];t2[[j]]=t2[[j]]+p1*p2[[k]],{k,0,2*198}],{j,1,steps}];
Clear[temp2];
Do[nup = lengthsup[[j]];ndown=lengthsdown[[j]];
If [nup>0,{sumval,temp2} =getIntegrate2[callup[[1;;2*198+1,j]],1,1,2*198+1];callup1[[j]]=sumval/t1[[j]]/dis[[j]],callup1[[j]] = 0];
If [ndown>0,{sumval,temp2} = getIntegrate2[calldown[[1;;2*198+1,j]],1,1,2*198+1];calldown1[[j]]=sumval/t2[[j]]/dis[[j]],calldown1[[j]]=0],{j,1,steps}];
Clear[temp2];
{confcall,temp2} =getIntegrate2[callup1,1,1,steps];{confput,temp2} =getIntegrate2[calldown1,1,1,steps];
confcall = confcall/steps;confput=confput/steps;Clear[temp2,callup,calldown,calldown1,callup1,t1,t2,valsupord,valsdownord,valsup,valsdown];
Do[valuesinupd[[j]]=newpassedthresholdsinup[[j]]/dis[[j]];
valuesindownd[[j]]=newpassedthresholdsindown[[j]]/dis[[j]];
nup
lengthsup[[j]];
valueforwardup[[j]]Re[Total[pdB[[indxup[[1;;nup,j]],j]]*newpassedthresholdforwardup[[indxup[[1;;nup,j]],j]]]/sumprob[[j]]];
ndown
lengthsdown[[j]];
valueforwarddown[[j]]Re[Total[pdB[[indxdown[[1;;ndown,j]],j]]*newpassedthresholdforwarddown[[indxdown[[1;;ndown,j]],j]]]/sumprob[[j]]],{j,1,steps}];
If[debugflag
1,Print["FinishedCalcFutureValues"]];
{dis,callup2,calldown2,confcall,confput,confcall10,confcall50,confcall60,confcall90,confput10,confput50,confput60,confput90,valueforwardup,valueforwarddown,valuesinupd,valuesindownd,temperatures,meanpath,futurearma1}]}
Compile->{Null}
googlePrices
(
All Rights Reserved - CopyRight 2012, Richard John Belshaw, AGU,AAAS,ASA,IEEE,ACS,ACM*)
ClearAll[googlePrices]; (* Google calls prices from last to first eg {last price,....,earliest price} )
googlePrices[stockTicker_String,startDate:(List|_String),endDate:(_List|_String),dataType_String:"Close"]:=Module[{price,volume,exportFormat"csv",newDataType},
( Importing data from Google
If[And @@ NumericQ /@ startDate,startDate=DateString[startDate,{"MonthNameShort"," ","Day",",","Year"}],startDate=DateString[startDate,{"MonthNameShort"," ","Day",",","Year"}]];
If[And @@ NumericQ /@ endDate,endDate=DateString[endDate,{"MonthNameShort"," ","Day",",","Year"}],endDate=DateString[endDate,{"MonthNameShort"," ","Day",",","Year"}]];
ds=ImportStringJoin[", "Data"]; )
newDataType=Switch[dataType,"Open"|"open","Open","High"|"high"|"hi"|"Hi","High","Low"|"low"|"Lo"|"lo","Low","Close"|"close","Close",
,dataType];
( ds=FinancialData[stockTicker,"OHLCV",{startDate,endDate}];
{ds, data}FinancialData[stockTicker,#,{startDate,endDate}]& /@ {"OHLCV",newDataType};;)
price=FinancialData[stockTicker,newDataType,{startDate,endDate}];
volume = FinancialData[stockTicker,"Volume",{startDate,endDate}]; {price,volume}]
{dddData,Volume}=googlePrices["DDD",{2011,1,1},{2013,2,24},"Close"]
DateListPlot[dddData]
GraphicsBox[{{}, {{},
{RGBColor[0.24720000000000014`, 0.24, 0.6], PointBox[CompressedData["
1:eJxVmHl0VOUZxi90O8Xtori0tTi4IO4DISRkgRvIRtZJwpaF5CaBAGHJDVkI
CcuFJJCEAMOaBAVH27qL49ZTNxzErR5tR0BRa3sGF8SAdjytG2rpzfO9D+fc
+Yfz49ue7/2e9/2+mzFVdYULhmua9swwTRv61/kZHd7BlMwTraMObBpvOGwd
cLht1IGME2fBgfccHmrOGDNhiMPDxw+m/NFKHnH8DzFDrN3mcP0QBsGhdxwe
QusncPScw4e/Hle2bTjYc/OEwZRvHHzpa8znm+nwmYK/jO5cAo4ccfitN53f
hWD9pwmcD3qMG2MGU6DcVnoLYjg+gfaHL5paPpTaPeejRF9UbB5/cTBlL0D
zu836O/Pc9jZbeuoL9EeaplIfeDovQ53nnYEvgK2x8S65gtmO+yoc/5H6W+K
lfUfUPoDsaI/Gfs33nDY2a0zo4pn5iQZfyVYa5jkio93n8OFQwNGYrz5Gtuf
Bvu/miTnlzER8dwbJ/u5PRZ6X47jfsHBL+NkvoJJ0HtVvKw/B2z1OewsPq7s
EnDgUDzjDQ6fjpf9rImD3ismy3rLwV6DfAAcOjhZ9tsCjp6aLHoGJsMPlyVw
PrAvOYH7S4D+hQ5fMxSuf4H1kwkSfy/Y0BOlfxPGWwmJ9IPSPz9R/NSr9G9N
lPNZhHiYFyWJ/04jfv44h3FcI8ChyiSuB45uTpL5N+B8PH9OknxJUH6YmMz4
34H4lieLvli0693JPH+0G086XIOf8vM/kyXemvJD2RSeH87bu3GKy59mcIpr
vP+DKZKvJzB/6OdTmV+3I77tDpfhgMH2o2w/7YX+41NF33fgyDBD/J8G1m91
OAY/VR8eMlz5FD5miJ/6lP5zhsx/i/LzTSmyP0PpLzrPmN/zYwr1QL9v7DT6
QeWfbxrnU/FunUY9GG98P437A1vXTWc8wIHc6S794ZXTJX73gbV7HR7q3boL
f2eVMmX7eBQVirjD33RxlTxwzTE03N3qvjFuhX6/5oq/t1/G+KZkcZ8wnh9
RZr4cZzSf1ea5N829LdeZX8NHIimsT7fgnjWpbMdbA6ky/4/RH4XTp36r8
8EW6K17RKzNk/yeVH/ZkSP9nlR9CGaxPiE9kMMMVL/3yTInXMaV/aqbL7+EX
MmX/5WjXTmXSL6o+XzpD9puPdjNpBvWr+NfMYP1Af8+n59vV/XFJFvUrf0zO
4vkjX4LVWbL+aOXnLVmu/VsXZrv8HZiUzfqFeITNbNd4bfN5Vvqfzmb9VfdH
TA7rg4rvvBypJ28of3flcP/K30/ksB6j3f4wx+VXvTTX5UejM1fa71H6H8uV
8W8r/e/nuuIb/lme63zMDXmsB1jf/0ie6OvBfkPv5tGf8I/nlnzWV/g5+GA+
44+kaP5rJdg/X/5cl7xqr6N81Ef2q1CN2tv++Q9cR36e3/wueJj3lDA+ZXe
AKXX0OrClzn6fuuQM6zE/3tawvJ6B/MKeT62E+kuZDrY7/6PYW8/1S+XVPE
fEf/8Iwi1nelv7GI9xn6e/cXMT/B5utF1KviWT/T5S/fnTO5X7Tbr8zkeBXv
f8+U++tNsNU/i/cHxgdemsX+yq9nZjFflV+vnM37UOXT7tmsR6oevDib8VV+
/Xy2+GWUyrdRcxgPsG/KHJe/9M/muPobI+fKeTSo+yxxrsS3Wt0XC+byvlL3
wSdzWW/A5sXFrveUP75Yzu8X6n1ZVczzUu/L3mJZfybYvqCE94PK/9gS5is4
UlHiej/pPSXi56txvxtPlfC9i/dBeEIp6wlYm1cq8z+M/t5NpcxvsPl4Ke8L
9Z74RynvB7CnpIznqd5rHWXcP9g+UCbrV4KD75UxP7F+ZPg8vm/RHnh4Huu
0vvOPL5n0a5p5ayv6n12c7nEf2w84vlAOfMDHD1S7nq/eX5i/xq8n3w3VjB/
wHZBBfWB9XAF8xFsnK1gvKHPut7kew4cyDPpJ/U+azGZL+p99q1Jv4P9YypZ
f8Ch7EqXnmhTJd9rSn+A7SPBwdFVzF9wJLOK9Vbpb6iif9Fu7KuiP8DWa1WM
v3r/plcz/oif16rm/Y12c2814wH2v1zNeKv38ZfVjB/Yt2w+30OYz+6bz/XA
wUPz6R9w5PR8+h8c2LWA7y/MFz64gO9z7E/7fAHjrfxwWQ3zE+PN5BreD8oP
z9XwfsN723OyhvUuEXr1hXxfg+2EhewPDs5fKPv/D9j4eCHf+2DrokW835Og
P24Rv5fA4cpF/B5Mhv7eRcwnsH/EYs4HDk1czPydCv3li1nfpkB/92L6H/19
Ty7m/Yr1IuNruT7a9bJa13rGxlrWE/S3grX83sN+Ah/U0j9gb/ESvpfAZvsS
xhfsf3QJ8w3xDR1fwvzA90x02FKZ70mwbS9lfVDxfWgpz0/pP0b+WOk/t5Tf
O9i/cdMyvh/A2pFl9BPi5f1xGes52By7nPmP/n4f+STW8/x9OfNRxfP75Txv
9Levq2P9NKA3t471GfNHVtYxX8HWN3X0Dzjgseh/cDjLYn4pPzRZ0v816PHe
bcl6dyk/L6e/oSeaEY944X5PCvqXfv33VXP+ZX+V+tZ79Gup62Q/AyBjboV
Es/12J81sELq+W9ToP8w2+9Be/iLFaxHaDeXNtAfYP+eBtnv6+DoYIPUz8un
Qe/ljYwv2oM7GyUfnwFHXmjk9xvW0081ip/Wgo1Lm+hvpTepie9VsPZck+zP
g/m8nzaxXii9lzTTX0rv5GbxSyM4VN3M9yLY91Gz5N/NYPvClbLfdKV/0kqp
hzco/eZKvl/A+uaVjAf2H/h1C/VPRzxjWlh/wVp5C88f/b1dLaK3HWw+wf4G
+ke9q0TPL8Ge0lX8fgP7Olfx/MH2Y6tkfAQcfJ+spSK+c1tdbG1olfh+hP6B
R9zt4XdbJb7FYG1YG/MV7F/Xxr+vgUMPtolfb02D/qNtct5XgT3/a+N7Ff19
41bz/gFH7lst/fvB+tur+fc/jDd+WM36mA79N6zh9xzaA/lrxN8fgM3v1vD7
Auy/di2/H8ChnLXMf6W3ea2sPxfr21+f76/iec061mNwZMY6id+1YL1xnew
B2zsXyf+Pqj8cLXNeqn8kGFLvn2l/FBvi/4S5Yc7bcZD+fkVW+LfBfakrnfl
o2/5elf+2v3rxV/TVX17ab2sfxQcObNe9C5V+bZkA99n4MDuDa56EX5xg/jt
AujTBjfI+fuU/lHt3I/Ktx3tjI+qF8+38++lmN/zWbvoKVL6R3ZI+z6lP7GD
+Yj59Wc7JB8/BRufdNBPiKd1cSf/nqD8HN/J/aJ/uKqT7yUV3xOdvJ/B/gs2
8v4Eh2I3cn6lv2Ij4wV/eHrYPgv+0edtYv1Vft20ifmp/Pr4Jn7Po90pALwv
wd6SLr4nlX87ujhe+fdAF7+XMqDvvS7mK9g3p5v1GP3t9d38XoHe4MPd9IvK
t3e6XfHUtR5+z6j4re3h39NVPXigh/eFqgdHeyQfngKbN27m+wAc/dNmnh/0
ecJs3494+M5u5ns3E3qv75Xz2on+wbxe5gP6G2/1SjyysT/r217Z37cq/8ds
ET+GweHsLYw3xmvNW2S9v4H9/93CfFHxHL2V711wNHMr30NKfwPbR6h479sq
9ekzzBf53TZ+2E/evo2vt9nQL+1jfECW3vZfggceHkb/Qb2TvfzPgX7+/zi
vzPg0CFyDjh62k/9YLt2u7Q/kgV9B7fzfsyGvs+3Sz6szoG+y3Zwv7mI3/Yd
9Hce4vf8DuYr2r0nyTFoN/WdEo/70e5P2Mn7AO2eZ3byPMG+j3fyvsuH3ot2
8XzBwbhd9AM4UrmL97kP8Rqxm+/ZAuiduFv8dgVYq9gt8agFe7t3y/hUcOhX
e5iv0BsdT9awnqdsj2t/vo17GH/Eyw7uob8wXr+jj+9NsFHcJ+tnYT6rvU/u
0/eV/kf7ZP52cPg4+z8LNmf3u/r77X7Wf8wXeqhf/L8Q7dFj/aLfVPrP9ct5
H4aeyP0DfC9hP/qRAfpL6f1xQPKhHv6wxu6V/qtyjP8DMOb/iA
"]]}, {}}, {}},
AspectRatio->0.6180339887498948,
AxesOrigin->{3515097600, 14.},
Frame->True,
FrameLabel->{None, None},
FrameTicks->{{Automatic, Automatic}, {{{3502828800, FormBox[""Jan"", TraditionalForm]}, {3518467200, FormBox[""Jul"", TraditionalForm]}, {3534364800, FormBox[""Jan"", TraditionalForm]}, {3550089600, FormBox[""Jul"", TraditionalForm]}, {3565987200, FormBox[""Jan"", TraditionalForm]}, {3581625600, FormBox[""Jul"", TraditionalForm]}, {3507926400, FormBox["""", TraditionalForm]}, {3513196800, FormBox["""", TraditionalForm]}, {3523824000, FormBox["""", TraditionalForm]}, {3529094400, FormBox["""", TraditionalForm]}, {3539548800, FormBox["""", TraditionalForm]}, {3544819200, FormBox["""", TraditionalForm]}, {3555446400, FormBox["""", TraditionalForm]}, {3560716800, FormBox["""", TraditionalForm]}, {3571084800, FormBox["""", TraditionalForm]}}, {{3502828800, FormBox["""", TraditionalForm]}, {3518467200, FormBox["""", TraditionalForm]}, {3534364800, FormBox["""", TraditionalForm]}, {3550089600, FormBox["""", TraditionalForm]}, {3565987200, FormBox["""", TraditionalForm]}, {3581625600, FormBox["""", TraditionalForm]}, {3507926400, FormBox["""", TraditionalForm]}, {3513196800, FormBox["""", TraditionalForm]}, {3523824000, FormBox["""", TraditionalForm]}, {3529094400, FormBox["""", TraditionalForm]}, {3539548800, FormBox["""", TraditionalForm]}, {3544819200, FormBox["""", TraditionalForm]}, {3555446400, FormBox["""", TraditionalForm]}, {3560716800, FormBox["""", TraditionalForm]}, {3571084800, FormBox["""", TraditionalForm]}}}},
GridLines->{{{3502828800, GrayLevel[0.8]}, {3518467200, GrayLevel[0.8]}, {3534364800, GrayLevel[0.8]}, {3550089600, GrayLevel[0.8]}, {3565987200, GrayLevel[0.8]}, {3581625600, GrayLevel[0.8]}}, None},
Method->{},
PlotRange->{{3515097600, 3.57048
^9}, {13.5, 69.8}},
PlotRangeClipping->True,
PlotRangePadding->{{1.107648*^6, 1.107648*^6}, {1.126, 1.126}},
Ticks->None]
optionPrices
(* All Rights Reserved - CopyRight 2012, Richard John Belshaw, AGU,AAAS,ASA,IEEE,ACS,ACM*)
ClearAll[optionPrices];
Compile -> {
optionPrices[tick1_String:"DDD",tick2_String:"DDD",dataType_String:"Close",StartDate:(List|_String),ExpDate:(_List|_String),strike:65,inflation_:0.02,interestRate_:0.02,dividend_:0.0,ppast_:0.5,paths_:1000,IndicatorPaths_Integer:1000,extendflag_Integer:1,forecastindicatorflag_Integer:0,currencyflag_Integer,debugflag_Integer:0];(0< ppast<= 1&&(extendflag=0||extendflag
1)&&(forecastindicatorflag
1 ||forecastindicatorflag
0)&&(currencyflag
1||currencyflag
0)&&(debugflag
1||debugflag
0)):=Module[{enddate,datatype,ds,TempData,expdate,s0,BDAYS,days,year,T,k1,infl,r,d,futuresteps,xa, xb, ratios, phaselist, extenlst, intlist, pdflist, outlst, difflist,xa3, xb3,ya3,yb3, ratios3, extenlst3, intlist3, pdflist3, outlst3, count1,count2,sumup,sumdown,rat1,rat2,pr,percChange,sumv1,intl1,tempd,sumv,intl,netGain,plotPrices,plotpdf,dB, pdB,i,j,count,dis,confcall10,confcall90,confput10,confput90, Tvalforwardup, Tvalforwarddown,valueforwardup,valueforwarddown, temperatures,meanpath,val1,val2,val3,val4,y,m,futurearma,futurearma1,conf01,conf10,conf25,conf50,conf75,stockTicker,conf90,conf99,temp,temp2,temp3,s1,s2,ratio33,ratio50, ratio80,ratio99,vol1,Volume,voldiff,datadiff,pi,indicator,plotFutureCall,plotIndicators,plotFuturePut,plotConfLimits,plotMeanPath,valuesinupd,valuesindownd,coeffs,volatility,vol2,plotVolatility,cubPCR1,data,plotFutureARMA,N1,sv,val0,val5,indicator2,temper2,neg50flag,neg80flag,neg98flag,numpts=8192,sumsame,count3,pr2,sma1,sma3,sma10,sma18,sma18v,ddd,ya,yb,temp4,temp5,temp6,temp7,temp8,temp9,temp10,temp11,plotIndicators2,sma1v,sma3v,sma10v,dB3,pdB3,meanpath3,meanpath10,meanpath18,meanpath3v,meanpath10v,meanpath18v,temperatures3,indicatorspredict,indicatorspredict2,numbersteps,numpts2,data2,meanpath2,meanpath22,indicator3,indicator4,meanpath1,meanpath1v,phasedistpdf3,phasedistpdf,phaseplot,confema310,confema390,confema3v10,confema3v90,delta1,delta2,conf60,confput60,indx1,indx2,confcall60,confcall50,confput50,val6,val7,val8,val9,indxfirst,indxlast,maxval,indicatorpredictplain,indicatorpredict2plain,forecastindicatordays,ratiolist,indxlist,MaxRatio,negFlag,ratio,confcall,confput,VolumePlot,PercentWinPlot,callup,calldown,CallupPlot,CalldownPlot,img,k,kmax},
(this code simply sets all the input text fields to blank upon pressing the calculateplots button *)
If [tick2 =
"",stockTicker = tick1, stockTicker = tick2];
(*this code just gets todays date in a particular format mm/dd/yyyy)
enddate=DateString[{"MonthName", "","Day","","Year"}];
(this code gets whether we are working with price data or volume data ie.datatype is either'hi','open','lo','close' or'vol')
datatype=Switch[dataType,"Open"|"open","Open","High"|"high"|"hi"|"Hi","High","Low"|"low"|"Lo"|"lo","Low","Close"|"close","Close",,dataType];
(this code sets the number of paths (a) to forecast to be at least 32768 if volume time series is chosen)(and the number of paths (a) to calculate is set to 16384)
(get the stock ticker symbol from the text input fieldd)
(get one years worth of data,either volume or price data)
(get the expiration date,or final date for forecasting to for the option or the forecast)
expdate=ExpDate;
(* get the last days price,or volume value )
( put it into yesterdays price or volume input field for user viewing, priceVolumeInputField=ToString[s0];)
( a special matlab function that computes the business days between todays date and the expiration date )
BDAYS=Cases[DayOfWeek @ Table[DatePlus[Take[DateList[],3],i],{i,0,DaysBetween[Take[DateList[],3],Take[DateList[expdate],3]]}],Except[Saturday|Sunday]];
(*compute the length of this array to get the total number of such days)
days=Length[BDAYS];
(get the number of business days in the current year)
year=Length[Cases[DayOfWeek /@ Table[DatePlus[{2013,1,1},i],{i,0,DaysBetween[{2013,1,1},{2013,12,31}]}],Except[Saturday|Sunday]]]-9;
(ratio of number of business days till expiration to total number of business days in current year)
T=days/year;
(get the strike price from the input text field)
k1=strike; (* 35.5 )
(*get the inflation rate from the input text field)
infl=inflation; (* 0.05 )
(*get the risk free interest rate from the input text field)
r=interestRate; (* 0.03 )
(*get the continuous dividend rate from the input text field)
d=dividend; (* 0.0 )
(*set the number of futuresteps to calculate to be the number of business days until expiration date)
futuresteps=days;
forecastindicatordays = futuresteps;
{TempData,vol1}googlePrices[stockTicker,StartDate,enddate,datatype];
data =Transpose[TempData][[2]];
s0=data[[-1]];
If [strike
0,k1=s0];
If[ currencyflag =
0,
Volume = Transpose[vol1][[2]];
voldiff = ConstantArray[0,Length[Volume]-1]];
datadiff = ConstantArray[0,Length[data]-1];
Do[datadiff[[i]] = data[[i+1]]-data[[i]],{i,1,Length[data]-1}];
If[currencyflag =0,
Do[voldiff[[i]] = Volume[[i+1]]-Volume[[i]],{i,1,Length[Volume]-1}]];
pi = ArcTan[1,1]4;
indicator =N[datadiff];
If[debugflag
1,Print["*1"]];
If [currencyflag =
0,
indicator2 N[voldiff/Max[Abs[voldiff]]*20.],
indicator2 = RandomReal[{-.01,.01},Length[data]]];
If[debugflag
1,Print["*2"]];
If [forecastindicatorflag
1,
numbersteps = forecastindicatordays,numbersteps=0];
sma1 = ConstantArray[0,Length[data]+numbersteps];
sma3 = ConstantArray[0,Length[data]+numbersteps];
sma10 = ConstantArray[0,Length[data]+numbersteps];
sma18 = ConstantArray[0,Length[data]+numbersteps];
sma1v = ConstantArray[0,Length[data]+numbersteps];
sma3v = ConstantArray[0,Length[data]+numbersteps];
sma10v = ConstantArray[0,Length[data]+numbersteps];
sma18v = ConstantArray[0,Length[data]+numbersteps];
Do [ sma1[[i+1]]
N[indicator[[i]]],{i,1,Length[indicator]}];
Do [ sma3[[i+3]]= N[1/2*Total[indicator[[i;;i+2]]]/3 + 1/2*sma3[[i-1+3]]],{i,1,Length[indicator]-2}];
Do[sma10[[i+9]] N[ 2/10*Total[indicator[[i;;i+8]]]/9 + 8/10*sma10[[i-1+9]]],{i,1,Length[indicator]-8}];
Do[sma18[[i+18]] = N[2/19*Total[indicator[[i;;i+17]]]/18+17/19*sma18[[i-1+18]]],{i,1,Length[indicator]-17}];
Do [ sma1v[[i+1]]
N[indicator2[[i]]],{i,1,Length[indicator2]}];
Do [ sma3v[[i+3]]= N[1/2*Total[indicator2[[i;;i+2]]]/3 + 1/2*sma3v[[i-1+3]]],{i,1,Length[indicator2]-2}];
Do[sma10v[[i+9]] N[ 2/10*Total[indicator2[[i;;i+8]]]/9 + 8/10*sma10v[[i-1+9]]],{i,1,Length[indicator2]-8}];
Do[sma18v[[i+18]] = N[2/19*Total[indicator2[[i;;i+17]]]/18+17/19*sma18v[[i-1+18]]],{i,1,Length[indicator2]-17}];
If [forecastindicatorflag
1,
meanpath1 = ConstantArray[0,numbersteps];
meanpath1v = ConstantArray[0,numbersteps];
meanpath2 = ConstantArray[0,numbersteps];
meanpath22 = ConstantArray[0,numbersteps];
meanpath3 = ConstantArray[0,numbersteps];
meanpath10 = ConstantArray[0,numbersteps];
meanpath18
ConstantArray[0,numbersteps];
meanpath3v = ConstantArray[0,numbersteps];
meanpath10v = ConstantArray[0,numbersteps];
meanpath18v = ConstantArray[0,numbersteps];
numpts2 Max[4*Length[indicator],8192];
indicator3 = sma3[[4;;Length[data]]];
indicator4 = sma3v[[4;;Length[data]]];
delta1 =10*Max[Abs[indicator3]];
If [delta1 =
0,delta1 = 500];
indicator3 = indicator3 + delta1;
delta2 = 10*Max[Abs[indicator4]];
If [delta2 = 0,delta2 = 500];
indicator4 = indicator4+delta2;
If[debugflag
1,Print["
*1"]];
{xa3, xb3, ratios3, extenlst3, intlist3, pdflist3, outlst3} = getData[indicator3,numpts2,extendflag,debugflag];
If[debugflag
1,Print["
*2"]];
{dB3, pdB3,phasedistpdf3} = generateMagPhase[ratios3, intlist3,pdflist3,outlst3, xa3, xb3, IndicatorPaths, numbersteps,1.0,debugflag];
If[debugflag
1,Print["
*3"]];
temperatures3 = generateTemperatures[dB3,indicator3[[-1]],IndicatorPaths,numbersteps];
If[debugflag
1,Print["
*4"]];
Do[ meanpath2[[j]] = Total[pdB3[[1;;IndicatorPaths,j]]*temperatures3[[1;;IndicatorPaths,j]]]/Total[pdB3[[1;;IndicatorPaths,j]]],{j,1,numbersteps}];
negFlag=0;If [Mean[meanpath2]
1,Print["
*5"]]; {xa3, xb3, ratios3, extenlst3, intlist3, pdflist3, outlst3} = getData[indicator4,numpts2,extendflag,debugflag]; dB3 = ConstantArray[0,{IndicatorPaths,numbersteps}]; pdB3 = ConstantArray[0,{IndicatorPaths,numbersteps}]; temperatures3 = ConstantArray[0,{IndicatorPaths,numbersteps}]; phasedistpdf3 = ConstantArray[0,Length[ratios]]; {dB3, pdB3,phasedistpdf3} = generateMagPhase[ratios3, intlist3,pdflist3,outlst3, xa3, xb3, IndicatorPaths, numbersteps,1.0,debugflag]; temperatures3 = generateTemperatures[dB3,indicator4[[-1]],IndicatorPaths,numbersteps]; Do[ meanpath22[[j]] = Total[pdB3[[1;;IndicatorPaths,j]]*temperatures3[[1;;IndicatorPaths,j]]]/Total[pdB3[[1;;IndicatorPaths,j]]],{j,1,numbersteps}]; negFlag=0;If[Mean[meanpath22]
1,Print["
*6"]]; confema310 = confema310 - delta1; confema390 = confema390-delta1; confema3v10 = confema3v10 - delta2; confema3v90 = confema3v90 - delta2; meanpath2 = meanpath2-delta1; meanpath22 = meanpath22-delta2; indicator3 = sma3[[1;;Length[data]]]; indicator3 =Join[indicator3,meanpath2]; indicator4 = sma3v[[1;;Length[data]]]; indicator4 = Join[indicator4,meanpath22]; sma3 = indicator3; Do[sma10[[i+2]] =N[1/2*Total[indicator3[[i;;i+2]]]/3 + 1/2*sma10[[i-1+2]]],{i,Length[data]-1,Length[indicator3]-2}]; Do[sma18[[i+5]] = N[2/7*Total[indicator3[[i;;i+5]]]/6+5/7*sma18[[i-1+5]]],{i,Length[data]-16,Length[indicator3]-5}]; sma3v = indicator4; Do[sma10v[[i+2]] =N[1/2*Total[indicator4[[i;;i+2]]]/3 + 1/2*sma10v[[i-1+2]]],{i,Length[data]-1,Length[indicator4]-2}]; Do[sma18v[[i+5]] = N[2/7*Total[indicator4[[i;;i+5]]]/6+5/7*sma18v[[i-1+5]]],{i,Length[data]-16,Length[indicator4]-5}]; meanpath3 = sma3[[Length[data]+1;;Length[data]+numbersteps]]; meanpath10=sma10[[Length[data]+1;;Length[data]+numbersteps]]; meanpath18 = sma18[[Length[data]+1;;Length[data]+numbersteps]]; meanpath3v = sma3v[[Length[data]+1;;Length[data]+numbersteps]]; meanpath10v = sma10v[[Length[data]+1;;Length[data]+numbersteps]]; meanpath18v = sma18v[[Length[data]+1;;Length[data]+numbersteps]]; ]; sma3 = sma3[[1;;Length[data]]]/Max[Abs[sma3[[1;;Length[data]]]]]*20; sma10 = sma10[[1;;Length[data]]]/Max[Abs[sma10[[1;;Length[data]]]]]*20; sma18 = sma18[[1;;Length[data]]]/Max[Abs[sma18[[1;;Length[data]]]]]*20; sma18v = sma18v[[1;;Length[data]]]/Max[Abs[sma18v[[1;;Length[data]]]]]*20; sma3v = sma3v[[1;;Length[data]]]/Max[Abs[sma3v[[1;;Length[data]]]]]*20; sma10v = sma10v[[1;;Length[data]]]/Max[Abs[sma10v[[1;;Length[data]]]]]*20; (*read into paths variable number of paths to calculate into future for many worlds probabilistic analysis..... ie.monte carlo) (set the number of points in the pdf to be 8192,for interpolation, numpts=8192;) (do the interpolation,and get the pdf and x array and phaselist and difflist) {xa, xb, ratios, extenlst, intlist, pdflist, outlst} = getData[data, numpts,extendflag,debugflag]; (do a quick analysis of ups versus downs in price or volume data) count1=0;count2=0;sumup=1;sumdown=1;sumsame =1; count3=0; ddd = year; If [Length[ratios] 1,count1=count1+1;sumup=sumup*ratios[[-i]],If[ratios[[-i]]<1,count2=count2+1;sumdown=sumdown*ratios[[-i]]]],{i,1,ddd}];
(rat1 is the geometric mean of ratios of number of up days in past sequence)
If [count1 ==0, rat1 = 0,
rat1 = sumup^(1/count1)
N];
(rat2 is the geometric mean of ratios of number of down days in past sequence)
If [count2
0,rat2=0,
rat2 = sumdown^(1/count2)
N];
(pr is probability of an up day in sequence compared to total number of days in sorted ratiolist)
If [count1+count2+count3
0,pr=0;pr2=0,
pr = count1
(count1+count2+count3);
pr2 = count2/(count1+count2+count3)];
(output these values for user viewing and compute the percentage increase or decrease out of 100 based on price or volume data from start of data to end of data (TempData(1)))
percChange=100*data[[-1]]data[[-ddd]]-100;
(set the focus to axes5) (erase previous data from axes)
(plot the price or volume data)
(plot the pdf,against its x axis values)
If [debugflag
1,Print["**1"]];
indx1 =matlabFind1[outlst,GreaterEqual, 0.9];
indx2 = matlabFind1[outlst,LessEqual,1.1];
indxfirst = indx1[[1]];
indxlast = indx2[[-1]];
If[debugflag
1,Print["
*1a"]];
{sumv1,intl1}=getIntegrate2[pdflist[[indxfirst;;indxlast]],1,1,indxlast-indxfirst+1];
If[debugflag
1,Print["
*1b"]];
tempd=pdflist[[indxfirst;;indxlast]]*outlst[[indxfirst;;indxlast]]/sumv1;
If[debugflag
1,Print["
*1c"]];
{sumv,intl}=getIntegrate2[tempd,1,1,indxlast-indxfirst+1];
If[debugflag
1,Print["
*1d"]];
Clear[indx1,indx2];
indx1 = matlabFind1[outlst,GreaterEqual,0.95];
If[debugflag
1,Print["
*1e"]];
indx2 = matlabFind1[outlst,LessEqual,1.05];
indxfirst = indx1[[1]];
indxlast = indx2[[-1]];
If[debugflag
1,Print["
*1f"]];
Clear[indx1];
indx1 = Reverse[Sort[pdflist[[indxfirst;;indxlast]]]];
maxval
indx1[[7]];
If[debugflag=1,Print["
*2"]];
(*compute the theoretical net gain and output it for user viewing,similarly output the expected value of the pdf)
If [rat1
0||rat2
0,netGain=percChange,
netGain
(rat1^(pr*year)rat2^(pr2*year)-1)*100];
(*plot the pdf,against its x axis values)
plotpdf=ListLinePlot[Transpose[{outlst[[indxfirst;;indxlast]],pdflist[[indxfirst;;indxlast]]}],Filling-> Axis,PlotRange-> {{.95,1.05},{0,maxval}},PlotLabel->"Probability Distribution",GridLines->Automatic,AxesLabel->{"price ratios","relative prob"}];
(clear the rest of the plots just for now)
(get the ratio probability ppast of number of phase points to calculate versus generic phase from the long distribution phaselistpdf2.)
(generate the magnitude and phase data,along with the probabilities of each ratio)
{dB, pdB,phasedistpdf} = generateMagPhase[ratios, intlist,pdflist,outlst, xa, xb, paths, futuresteps,ppast,debugflag];
{coeffs,volatility,vol2}getPercentErrorCubic[data,{{}}];
cubPCR1=ConstantArray[0,Length[data]];
Do [ cubPCR1[[j]] =coeffs[[1]]j^3+coeffs[[2]]*j^2+coeffs[[3]]*j+coeffs[[4]],{j,1,Length[data]}];
temp = ConstantArray[{},Length[data]];
temp4 = ConstantArray[{},Length[data]];
temp5
ConstantArray[{},Length[data]];
temp6 = ConstantArray[{},Length[data]];
temp7 = ConstantArray[{},Length[data]];
temp8 = ConstantArray[{},Length[data]];
temp9 = ConstantArray[{},Length[data]];
temp10 = ConstantArray[{},Length[data]];
temp11 = ConstantArray[{},Length[data]];
Do [temp[[j]] = {Transpose[TempData][[1,j]],cubPCR1[[j]]}; temp4[[j]] = {Transpose[TempData][[1,j]],sma3[[j]]};temp5[[j]] = {Transpose[TempData][[1,j]],sma10[[j]]};temp6[[j]]= {Transpose[TempData][[1,j]],sma18[[j]]};temp7[[j]]= {Transpose[TempData][[1,j]],sma1[[j]]};temp8[[j]]= {Transpose[TempData][[1,j]],sma18v[[j]]};temp9[[j]]= {Transpose[TempData][[1,j]],sma10v[[j]]};temp10[[j]]= {Transpose[TempData][[1,j]],sma3v[[j]]};temp11[[j]]= {Transpose[TempData][[1,j]],sma1v[[j]]},{j,1,Length[data]}];
plotPrices=Show[DateListPlot[TempData,Joined->True,PlotRange-> Automatic,GridLines->Automatic,AxesLabel->{"Past Days","underlying"},PlotLabel->"Past years data for "<>stockTicker,RotateLabel->True],DateListPlot[temp,Joined-> True,GridLines->Automatic,PlotStyle-> Directive[Red ]] ];
plotIndicators Show[
DateListPlot[temp6,PlotRange-> {-20,20},Joined-> True,PlotStyle-> Directive[Blue],PlotLabel-> "Price Delta Indicators 3(Magenta) 9(Green) 18(Blue) EMA day"],
DateListPlot[temp4,Joined-> True,PlotStyle->Directive[Magenta]],DateListPlot[temp5,Joined-> True,PlotStyle->Directive[Green]]
];
plotIndicators2 = Show[
DateListPlot[temp9,PlotRange-> {-20,20},Joined-> True,PlotStyle-> Directive[Green],PlotLabel-> "Volume Delta Indicators 3(Magenta) 9(Green) 18(Blue) EMA day"],
DateListPlot[temp10,Joined-> True,PlotStyle->Directive[Magenta]],
DateListPlot[temp8,Joined-> True,PlotStyle->Directive[Blue],PlotLabel-> "Volume/Price Indicators, 1,3,10,18 day"]];
If [forecastindicatorflag ==1,
temp2 = ConstantArray[{},numbersteps];
temp3 = ConstantArray[{},numbersteps];
temp5 = ConstantArray[{},numbersteps];
temp6 = ConstantArray[{},numbersteps];
temp7 = ConstantArray[{},numbersteps];
temp9 = ConstantArray[{},numbersteps];
temp10 = ConstantArray[{},numbersteps];
temp11 = ConstantArray[{},numbersteps];
Do [temp2[[j]] = {j,confema310[[j]]};temp3[[j]]
{j,confema390[[j]]};temp5[[j]] = {j,meanpath3[[j]]};temp6[[j]]= {j,meanpath10[[j]]};temp7[[j]]= {j,meanpath18[[j]]};temp9[[j]]= {j,meanpath3v[[j]]};temp10[[j]]= {j,meanpath10v[[j]]};temp11[[j]]= {j,meanpath18v[[j]]},{j,1,numbersteps}];
indicatorspredict = ListLinePlot[{temp5,temp2,temp3,temp6,temp7},Joined-> True,GridLines ->Automatic,PlotRange->{Min[Transpose[temp2][[2]]],Max[Transpose[temp3][[2]]]},PlotStyle->{Magenta,Darker[Cyan],Darker[Cyan],Green,Blue},PlotLabel-> "Price Delta Forecast,n 3(magenta),9(Green),18(Blue);50% conf lim Dark Cyan"];
indicatorpredictplain = ListLinePlot[{temp5,temp6,temp7},Joined-> True,GridLines ->Automatic,PlotRange-> Automatic,PlotStyle->{Magenta,Green,Blue},PlotLabel-> "Price Delta Forecast,n 3magenta,9Green,18Blue"];
indicatorpredict2plain = ListLinePlot[{temp9,temp10,temp11},Joined-> True,GridLines ->Automatic,PlotRange->Automatic,PlotStyle->{Magenta,Green,Blue},PlotLabel-> "Volume Delta Forecast,n 3magenta,9Green,18Blue"];
Do [temp2[[j]] = {j,confema3v10[[j]]};temp3[[j]]{j,confema3v90[[j]]},{j,1,numbersteps}];
indicatorspredict2 = ListLinePlot[{temp9,temp2,temp3,temp10,temp11},Joined-> True,GridLines ->Automatic,PlotRange->{Min[Transpose[temp2][[2]]],Max[Transpose[temp3][[2]]]},PlotStyle->{Magenta,Darker[Cyan],Darker[Cyan],Green,Blue},PlotLabel-> "Volume Delta Forecast,n 3(magenta),9(Green),18(Blue);50% conf lim Dark Cyan"];
];
(*Do[ temp8[[j]] = {j,phasedistpdf[[j]]},{j,1,Length[phasedistpdf]}];
phaseplot = ListPlot[temp8,PlotStyleDirective[Blue],PlotRange Automatic,PlotLabel{"Phase distribution by delta day from present"}];
*)
data2 = ConstantArray[0,Length[data]];
Do[ data2[[j]] = data[[j]] - cubPCR1[[j]],{j,1,Length[data]}];
N1 = Floor[Log[Length[data2]]/Log[2.0]]+1;
temp3 = Join[data2, ConstantArray[0,2^N1-Length[data2]]];
y=Fourier[temp3]/2^(N1-1);
futurearma=ConstantArray[0,futuresteps];(
Sinusoidal path )
Do[futurearma[[i]]=futurearma[[i]]+Re[Abs[y[[j]]]*Exp[I (i+Length[data2])/2^N1(j-1)+ArcTan[Re[y[[j]]],Im[y[[j]]]]]],{i,1,futuresteps},{j,2,2^(N1-1)}];
(calculate the option plot values and confidence limits and most probable path into the future)
{dis,callup,calldown,confcall,confput,confcall10,confcall50,confcall60,confcall90,confput10,confput50,confput60,confput90, valueforwardup, valueforwarddown,valuesinupd,valuesindownd, temperatures,meanpath,futurearma1} = calculateFutureValues[dB, pdB,futurearma, s0, k1, r, d, infl, T,debugflag];
plotFutureARMA = ListLinePlot[futurearma1,AxesLabel->{"days ahead","Fut sinusoidal"},PlotLabel->"Most probable pathn w
Sinusoidal volatility overlay",GridLines->Automatic,RotateLabel->True];
(calculate the limits and most probable values of the call and put options,output them for user viewing to text fields)
val0=0;val1 =0;val2=0;val3=0;val4=0;val5=0;
Do[val0=val0+valueforwardup[[j]]dis[[j]]; val1 = val1 + confcall10[[j]]/dis[[j]];val2 = val2 + confcall90[[j]]/dis[[j]];
val3 = val3 + confput10[[j]]/dis[[j]];val4 = val4+confput90[[j]]/dis[[j]];val5=val5+valueforwarddown[[j]]/dis[[j]],{j,1,futuresteps}];
val0 = val0/futuresteps;val5=val5/futuresteps;
val1 = val1/futuresteps;val2 = val2/futuresteps;
val3 = val3/futuresteps;val4=val4/futuresteps;
(calculate confidence levels for Risk analysis)
val6=0;val7=0;val8=0;val9=0;
Do[val6 = val6+confcall60[[j]]/dis[[j]];val7=val7+confput60[[j]]/dis[[j]];val8=val8+confcall50[[j]]/dis[[j]];val9=val9+confput50[[j]]/dis[[j]],{j,1,futuresteps}];
val6=val6/futuresteps;
val7=val7/futuresteps;
val8 = val8/futuresteps;
val9=val9/futuresteps;
If [debugflag
1,Print[Dimensions[temperatures]]];
negFlag=0;If [Mean[meanpath]{Blue, Red,Red},AxesLabel->{"days ahead","dollars"},GridLines->Automatic,PlotLabel->"Present Value: Call Option;n80% conf lim Red",RotateLabel->True];
(* plot future value of put option, plot confidence limits as well )
plotFuturePut=ListLinePlot[{valueforwarddown,confput10,confput90},PlotStyle->{Blue, Red,Red},GridLines->Automatic,AxesLabel->{"days from today's date","dollars"},PlotLabel->"Present Value: Put Option;n80% conf lim Red",RotateLabel->True];
( if volume data,plot log of confidence limits into future, otherwise plot confidence limits not in log form )
If[debugflag
= 1,Print["Startcubic"]];
{coeffs,volatility,vol2} = getPercentErrorCubic[data,temperatures];
If[debugflag=1,Print["Endcubic"]];
plotVolatility = ListLinePlot[volatility,AxesLabel->{"days ahead","Volatility"},GridLines->Automatic,PlotLabel->"Volatility Plot"];
plotConfLimits=ListLinePlot[{meanpath,conf50,conf01,conf99,conf25,conf75,conf10,conf90},PlotStyle->{Magenta,Blue, Darker[Green],Darker[Green],Darker[Cyan],Darker[Cyan],Red,Red},AxesLabel->{"future days","underlying $"},PlotLabel->"Confidence Lines",GridLines->Automatic,PlotLegends->{"mean path","50%", "1%","99%","25%","75%","10%","90%"},RotateLabel->True,ImageSize->{300,270}];
(
plot the most probable path, if price data,plot the 80% confidence limits as well )
plotMeanPath=ListLinePlot[{meanpath,conf10,conf90},PlotStyle->{Blue, Red,Red},GridLines->Automatic,PlotLabel->"Most Probable Path;n 80% conf limits Red",AxesLabel->{"future days","prob Fut Val"},RotateLabel->True];
VolumePlot = DateListPlot[vol1,Joined-> False,PlotStyle-> Red,PlotLabel-> "Volume on a given date",GridLines-> Automatic];
Clear[temp3];
temp3 = ConstantArray[{},198];
Do[temp3[[j]] = {indxlist[[-j]],ratiolist[[-j]]},{j,1,198}];
PercentWinPlot = ListLinePlot[temp3,PlotStyle-> Red,PlotLabel-> "% Win Chance versus probability of confidence level",GridLines-> Automatic,PlotRange-> Automatic,Joined-> True];
If [debugflag
1,Print[Dimensions[callup]]];
N1=50;
img=ConstantArray[{},Floor[(2*198+1-5+1)/4]+1];
k=1;
Do[i=j;If[i>199,i=2*198+1-i];img[[k]]=ListLinePlot[callup[[j,1;;futuresteps]],PlotStyle-> RGBColor[1,1-N[(N1+i)
(198+1+N1)],1-N[(N1+j)(198+1+N1)]],GridLines-> Automatic,Joined-> True];k=k+1,{j,5,2*198+1,4}];
kmax k-1;
CallupPlot Show[ ListLinePlot[callup[[1,1;;futuresteps]],PlotStyle-> RGBColor[1,1-N[N1
(198+1+N1)],1-N[N1/(198+1+N1)]],PlotRange-> {Min[callup[[1;;2*198+1,1;;futuresteps]]],Max[callup[[1;;2*198+1,1;;futuresteps]]]},PlotLabel->"Call values across days for different conf levels",GridLines-> Automatic,Joined-> True],
Table[img[[k]],{k,kmax}],ListLinePlot[valueforwardup,PlotStyle-> Blue,Joined-> True]];
Clear[img];
img=ConstantArray[{},Floor[(2*198+1-5+1)4]+1];
k=1;Do[i=j;If[i>199,i=2*198+1-i];img[[k]]=ListLinePlot[calldown[[j,1;;futuresteps]],PlotStyle -> RGBColor[1-N[(N1+i)(198+1+N1)],1-N[(N1+i)(198+1+N1)],1],GridLines-> Automatic,Joined-> True];k=k+1,{j,5,2*198+1,4}];
kmax=k-1;
CalldownPlot=Show[ ListLinePlot[calldown[[1,1;;futuresteps]],PlotRange-> {Min[calldown[[1;;2*198+1,1;;futuresteps]]],Max[calldown[[1;;2*198+1,1;;futuresteps]]]},PlotLabel->"Put values across days for different conf levels",
PlotStyle -> RGBColor[1-N[N1(198+1+N1)],1-N[N1/(198+1+N1)],1],GridLines-> Automatic,Joined-> True],Table[img[[k]],{k,kmax}],ListLinePlot[valueforwarddown,PlotStyle-> Red,Joined-> True]];
Clear[img];
Labeled[Framed[Column[{If [forecastindicatorflag
=1,GraphicsGrid[{{plotPrices,plotpdf},{plotIndicators,plotFutureARMA},{plotIndicators2,plotMeanPath},{plotFutureCall,indicatorpredictplain},{plotFuturePut,indicatorpredict2plain},{plotVolatility, indicatorspredict},{plotConfLimits,indicatorspredict2},{VolumePlot,PercentWinPlot},{CallupPlot,CalldownPlot}},ImageSize-> {800,2200}],
GraphicsGrid[{{plotPrices,plotpdf},{plotIndicators,plotFutureARMA},{plotIndicators2,plotMeanPath},{plotFutureCall,plotFuturePut},{plotVolatility,plotConfLimits},{VolumePlot,PercentWinPlot},{CallupPlot,CalldownPlot}},ImageSize-> {800,1800}]],
Grid[{{Labeled[Framed[ToString[stockTicker]],"Stock Ticker Symbol", Top],Labeled[Framed[ToString[dataType]],"Data Type", Top],Labeled[Framed[ToString[enddate]],"Today's date", Top],Labeled[Framed[ToString[ExpDate]],"Expiration Date", Top]},{Labeled[Framed[strike],"Strike Price", Top],Labeled[Framed[interestRate],"Risk Free Rate", Top],Labeled[Framed[inflation],"Inflation rate", Top],Labeled[Framed[dividend],"Dividend Yield", Top]},{Labeled[Framed[ppast],"ppast", Top],Labeled[Framed[paths],"Number of paths", Top],Labeled[Framed[percChange],"Percentage Increase", Top],Labeled[Framed[netGain],"Theoretical Gain Index", Top]},{Labeled[Framed[s0],"Yesterday's price", Top],Labeled[Framed[ratio],"% chance win ALL conf lev", Top]},{Labeled[Framed[ratio99],"% chance win 99% time", Top],Labeled[Framed[ratio80],"% chance win 80% time", Top],Labeled[Framed[ratio50],"% chance win 50% time", Top],Labeled[Framed[ratio33],"% chance win 33% time", Top]},{Labeled[Framed[N[pr]],"Prob of Ups", Top],Labeled[Framed[N[rat1-1]],"Mean Up", Top],Labeled[Framed[N[1-rat2]],"Mean Down", Top],Labeled[Framed[sumv],"mean value of Distbn", Top]},
{Labeled[Framed[val0],"Most Probable Call value", Top],Labeled[Framed[confcall],"Call Value, All conf lev",Top],Labeled[Framed[val5],"Most Probable Put value", Top],Labeled[Framed[confput],"Put Value, All conf lev",Top]},
{Labeled[Framed[val1],"80% lower limit of Call", Top],Labeled[Framed[val8],"50% median limit of Call",Top],Labeled[Framed[val6],"60% upper limit of Call",Top],Labeled[Framed[val2],"80% upper limit of Call", Top]},{Labeled[Framed[val3],"80% lower limit of Put", Top],Labeled[Framed[val9],"50% median limit of Put",Top],Labeled[Framed[val7],"60% upper limit of Put",Top],Labeled[Framed[val4],"80% upper limit of Put", Top]}
}]}]]
,Style["Option Prices, Future Values, Mean Path and Sinusoidal Path for "<>ToString[stockTicker]<>" from "<>ToString[enddate]<>" to "<>ToString[ExpDate],Brown,14], Top,Alignment->Center]]}
Compile->{Null}
Compile->{Null}
Compile->{Null}
Examples
( All Rights Reserved - CopyRight 2012, Richard John Belshaw, AGU,AAAS,ASA,IEEE,ACS,ACM*)
optionPrices[tick1_String:"DDD",tick2_String:"DDD",dataType_String:"Close",StartDate:(_List|_String),ExpDate:(_List|_String),strike:65,inflation_:0.02,interestRate_:0.02,dividend_:0.0,ppast_:0.5,paths_:1000,IndicatorPaths_Integer:1000,extendflag_Integer:1,forecastindicatorflag_Integer:0,currencyflag_Integer,debugflag_Integer:0]
optionPrices["AMGN","","Close","Jan 1,2012","July 1,2013",0,0.02,0.01,0.0,0.5,8192,3000,1,1,0,1]
WolframAlpha["AMGN",{{"FundamentalsAndFinancials",1},"ComputableData"}]
temp =Import["Symbols.csv"]
(Symbol
CRAY
OCN
HCI
BOFI
AZZ
GEOS
PRAA
AIRM
DHI
ACT
GIL
TMO
HOMB
WRLD
AMG
WAB
CPA
HCC
FCFS
UTHR
VRX
WDC
NDSN
COO
VPFG
)
optionPrices[tick1_String:"DDD",tick2_String:"DDD",dataType_String:"Close",StartDate:(List|_String),ExpDate:(_List|_String),strike:65,inflation_:0.02,interestRate_:0.02,dividend_:0.0,ppast_:0.5,paths_:1000,IndicatorPaths_Integer:1000,extendflag_Integer:1,forecastindicatorflag_Integer:0,currencyflag_Integer,debugflag_Integer:0]
Manipulate[optionPrices[tick1,tick2,dataType,StartDate,ExpDate,strike,inflation,interestRate,dividend,ppast,paths,indicatorPaths,extendflag,forecastindicatorflag,currencyflag,0],{{tick1,"CRAY","Stock Ticker"},{"CRAY","OCN","HCI","BOFI","AZZ","GEOS","PRAA","AIRM","DHI","ACT","GIL","TMO","HOMB","WRLD","AMG","WAB","CPA","HCC","FCFS","UTHR","VRX","WDC","NDSN","COO","VPFG","AMGN"},ControlType->PopupMenu },{{tick2,"","Input field for arbitrary Stock Ticker or currency eg. USD/EUR"},ControlType->InputField},
{{dataType,"Close","Data Type"},{"Open","High","Low","Close"},ControlType->PopupMenu},{{StartDate,"Jan 1,2012","Start Date for data"},ControlType->InputField},
{{ExpDate,"Dec 31,2013","Expiry Date"},ControlType->InputField},{{strike,0,"Strike"},ControlType->InputField},{{inflation,0.02,"Inflation"},ControlType->InputField},{{interestRate,0.02,"Interest Rate"},ControlType->InputField},{{dividend,0.0,"Dividend"},ControlType->InputField},{{ppast,0.5,"Past Prob"},ControlType->InputField},{{paths,8192,"No of Paths"},ControlType->InputField},{{indicatorPaths,3000,"No of Indicator Paths"},ControlType->InputField},{{extendflag,1,"Extend ratios by order statistics (0 or 1)"},{0,1},ControlType->Setter},
{{forecastindicatorflag,1,"Forecast indicators (0 or 1)"},{0,1},ControlType->Setter},{{currencyflag,0,"Currency (0 or 1)"},{0,1},ControlType->Setter},ControlPlacement->Top,SaveDefinitions->True,SynchronousUpdating->False]