Showing posts with label Black-Scholes formula-R. Show all posts
Showing posts with label Black-Scholes formula-R. Show all posts

Sunday, March 22, 2020

Black-Scholes formula-R

 Black-Scholes formula-R

> BlackScholes <- function(TypeFlag = c("c", "p"), S, X, Time, r, b, sigma) { TypeFlag = TypeFlag[1] d1 = (log(S/X) + (b + sigma * sigma/2) * Time)/(sigma * sqrt(Time)) d2 = d1 - sigma * sqrt(Time) if (TypeFlag == "c") price = S * exp((b - r) * Time) * pnorm(d1) - X * exp(-r * Time) * pnorm(d2) else if (TypeFlag == "p") price = X * exp(-r * Time) * pnorm(-d2) - S * exp((b - r) * Time) * pnorm(-d1) param <- list(TypeFlag = TypeFlag, S = S, X = X, Time = Time, r = r, b = b, sigma = sigma) ans <- list(parameters = param, price = price, option = "Black Scholes") class(ans) <- c("option", "list") ans }


> library("ggplot2", lib.loc="~/R/win-library/3.6")
> x <- seq(-4, 4, l = 50)
> y <- x
> myf <- function(x, y) {
+     sin(x) + cos(y)
+ }
> z <- outer(x, y, FUN = myf)
> persp(x, y, z, theta = 45, phi = 45, shade = 0.45)
> Cap = floor(Capitalization/1000)
> Cap
                2003  2004 2005  2006  2007 2008
Euronext US     1328 12707 3632 15421 15650 9208
TSX Group        888  1177 1482  1700  2186 1033
Australian SE    585   776  804  1095  1298  683
Bombay SE        278   386  553   818  1819  647
Hong Kong SE     714   861 1054  1714  2654 1328
NSE India        252   363  515   774  1660  600
Shanghai SE      360   314  286   917  3694 1425
Tokyo SE        2953  3557 4572  4614  4330 3115
BME Spanish SE   726   940  959  1322  1781  948
Deutsche Boerse 1079  1194 1221  1637  2105 1110
London SE       2460  2865 3058  3794  3851 1868
Euronext EU     2076  2441 2706  3712  4222 2101
SIX SE           727   826  935  1212  1271  857
> barplot(t(Cap)/1e+06, beside = TRUE, las = 2, ylab = "Capitalization [Mio USD]")
> title(main = "Major Stock Markets")
>
 function


> mtext(side = 3, "2003 - 2008")
> barplot(Cap/1e+06, beside = TRUE, ylab = "Capitalization [Mio USD]")
> palette(rainbow(13, s = 0.6, v = 0.75))
> stars(t(log(Cap)), draw.segments = TRUE, ncol = 3, nrow = 2,
+       key.loc = c(4.6, -0.5), mar = c(15, 0, 0, 0))
> mtext(side = 3, line = 2.2, text = "Growth and Decline of Major Stock Markets",
+       cex = 1.5, font = 2)
> abline(h = 0.9)
calldelta(100, 100, 0.30, (15/365), 0.02, 0.02)
[1] 0.5117085
> callgreek("delta", 100, 100, 0.20, (45/365), 0.02, 0.02)
[1] 0.5127391
> callgreek("delta", 100, 100, 0.20, (15/365), 0.02, 0.02)
[1] 0.5076694
> callgreek("gamma", 100, 100, 0.20, (45/365), 0.02, 0.02)
[1] 0.05663458
> callgreek("gamma", 100, 100, 0.20, (15/365), 0.02, 0.02)
[1] 0.09829574
> callpremium(100, 100, 0.20, (45/365), 0.02, 0.02)
[1] 2.794086
> callrho(100, 100, 0.20, (45/365), 0.02, 0.02)
[1] 0.05976964
> calltheta(100, 100, 0.20, (45/365), 0.02, 0.02)
[1] -0.02979083
> dv(s = 100, x1 = 90, x2 = 95, x3 = 105, x4 = 110, t = 0.08, r = 0.02, sigma = 0.2, vol = 0.3)
                    V1
Spot            100.00
Long.Put         90.00
Short.Put        95.00
Short.Call      105.00
Long.Call       110.00
Lower.BE         94.01
Higher.BE       105.99
Prob.BE           0.45
Prob.Max.Profit   0.38
Prob.Max.Loss     0.32
> iv.calc(type = "call", price = 42.93, s = 100, x = 100, t = (15/365), r = 0.02, d = 0)
[1] 5
> lambda(type = "put", s = 100, x = 100, sigma = 0.15, t = 45/365, r = 0.02)
[1] -23.80998
> opteval(100, 100, 0.20, (45/365), 0.02, 0.02)
  Option  Premium      Delta      Gamma
1   Call 2.794086  0.5127391 0.05663458
2    Put 2.794086 -0.4847982 0.05663458
       Vega       Theta         Rho
1 0.1396469 -0.02979083  0.05976964
2 0.1396469 -0.02979083 -0.06321441
> opteval(100, 9200, 0.20, (45/365), 0.02, 0.02)
  Option  Premium      Delta Gamma Vega
1   Call    0.000  0.0000000     0    0
2    Put 9077.589 -0.9975373     0    0
     Theta       Rho
1 0.000000   0.00000
2 4.034484 -11.31453
> opteval(100, 100, 0.20, (45/365), 0.02, 0.02)
  Option  Premium      Delta      Gamma
1   Call 2.794086  0.5127391 0.05663458
2    Put 2.794086 -0.4847982 0.05663458
       Vega       Theta         Rho
1 0.1396469 -0.02979083  0.05976964
2 0.1396469 -0.02979083 -0.06321441
> plotbearcall(s= 100, x1 = 95, x2 = 105, t = (45/365), r = 0.02,sigma = 0.20, sigma2 = 0.20, d = 0, ll = 0.75, ul = 1.25)
> plotbearput(s= 100, x1 = 95, x2 = 105, t = (45/365), r = 0.02,
+             sigma = 0.20, sigma2 = 0.20, d = 0, ll = 0.75, ul = 1.25)
> plotbullcall(s= 100, x1 = 95, x2 = 105, t = (45/365), r = 0.02,
+              sigma = 0.20, sigma2 = 0.20, d = 0, ll = 0.75, ul = 1.25)
> plotbullput(s= 100, x1 = 95, x2 = 105, t = (45/365), r = 0.02,
+             sigma = 0.20, sigma2 = 0.20, d = 0, ll = 0.75, ul = 1.25)
> plotdv(s= 100, x1 = 90, x2 = 95, x3 = 105, x4 = 110, t = (45/365), r = 0.02, sigma = 0.20)
> plotvertical("call", 100, 90, 110, (45/365), 0.02, 0.20)
> prob.above(spot = 100, lower = 110, mean = 0, dsd = 0.01, dte = 45)
[1] 0.06801856
> prob.above(spot = 100, mean = 0, dsd = 0.01, dte = 45, p = 0.75, quantile = TRUE)
  probability percent.change    price
1        0.75    -0.04524615 95.47539
> prob.below(spot = 100, upper = 110, mean = 0, dsd = 0.01, dte = 45)
[1] 0.9319814
> prob.below(spot = 100, mean = 0, dsd = 0.01, dte = 45, p = 0.75, quantile = TRUE)
  probability percent.change    price
1        0.75     0.04524615 104.5246
> prob.btwn(spot = 100, mean = 0, dsd = 0.01, dte = 45, p = 0.75, quantile = TRUE)
  probability percent.change     price
1        0.75    -0.07716778  92.28322
2        0.75     0.07716778 107.71678
> putdelta(100, 0.20, (45/365), 0.02, 0.02)
[1] 0
> puteval(100, 100, 0.20, (45/365), 0.02, 0.02)
   Premium      Delta      Gamma      Vega
1 2.794086 -0.4847982 0.05663458 0.1396469
        Theta         Rho
1 -0.02979083 -0.06321441
> putgreek("vega", 100, 100, 0.20, (45/365), 0.02, 0.02)
[1] 0.1396469
> putpremium(100, 100, 0.20, (45/365), 0.02, 0.02)
[1] 2.794086
> putrho(100, 100, 0.20, (45/365), 0.02, 0.02)
[1] -0.06321441
> puttheta(100, 100, 0.20, (45/365), 0.02, 0.02)
[1] -0.02979083
> r.cont(0.12, 2)
[1] 0.1165378
> vertical("call", s = 100, x1 = 90, x2 = 110, t = (45/365), r =  0.025, sigma = 0.20, vol = 0.25)
                    V1
Spot            100.00
Short.Strike     90.00
Long.Strike     110.00
Max.Profit       10.13
Max.Loss         -9.87
Breakeven       100.13
Prob.BE           0.50
Prob.Max.Profit   0.17
Prob.Max.Loss     0.17
Initial.DC       10.13
> AsianCall(T=1,d=12,K=100,r=0.05,sigma=0.2,S0=100,method=c("best"),
+           sampling=c("QMC"),metpar=list(maxiter=100,tol=1.e-14,cvmethod="splitting"),
+           sampar=list(nout=50,n=2039,a=1487,baker=TRUE,genmethod="pca"))
          result error estimate
price 6.15604075   3.336613e-08
delta 0.59385393   2.731083e-09
gamma 0.03057798   6.326157e-10
> AsianCall_AppLord(T = 1, d = 12, K = 100, r = 0.05, sigma = 0.2, S0 = 100)
[1] 6.15604
> AsianCall(T=1,d=12,K=170,r=0.05,sigma=0.2,S0=100,method="best",
+           sampling="MC",metpar=list(maxiter=100,tol=1.e-14,np=1000),
+           sampar=list(n=10^5))
            result error estimate
price 1.522912e-04   3.763444e-09
delta 5.091253e-05   8.851702e-10
gamma 1.578444e-05   1.759080e-10
> AsianCall_AppLord(T = 1, d = 12, K = 170, r = 0.05, sigma = 0.2, S0 = 100)
[1] 0.0001524371
> AsianCall_AppLord(T = 1, d = 12, K = 100, r = 0.05, sigma = 0.2, S0 = 100)
[1] 6.15604
> BS_EC(K=100, r = 0.05, sigma = 0.2, T = 0.25, S0 = 100)
     price      delta      gamma 
4.61499713 0.56946018 0.05694602 
> BS_EP(K=100, r = 0.05, sigma = 0.2, T = 0.25, S0 = 100)
      price       delta       gamma 
 3.37277718 -0.43053982  0.05694602 
> AsianCall(T=1,d=12,K=100,r=0.05,sigma=0.2,S0=100,method=c("best","naive"),
+           sampling=c("QMC","MC"),
+           metpar=list(maxiter=100,tol=1.e-14,cvmethod="splitting"),
+           sampar=list(nout=50,seq.type="korobov",n=2039,a=1487,
+                       baker=TRUE,genmethod="pca"))
          result error estimate
price 6.15604079   2.957191e-08
delta 0.59385393   2.440579e-09
gamma 0.03057798   5.350724e-10
> AsianCall(T=1,d=12,K=100,r=0.05,sigma=0.2,S0=100,method="best",
+           sampling="MC",metpar=list(maxiter=100,tol=1.e-14,np=1000),
+           sampar=list(n=10^5))
          result error estimate
price 6.15604075   1.945657e-07
delta 0.59385394   1.369256e-08
gamma 0.03057798   3.526989e-09
> AsianCall(T=1,d=12,K=100,r=0.05,sigma=0.2,S0=100,method="naive",
+           sampling="QMC",
+           sampar=list(nout=50,n=2039,a=1487,baker=TRUE,genmethod="pca"))
          result error estimate
price 6.15618443   0.0007920166
delta 0.59392376   0.0004117103
gamma 0.03059009   0.0001390430
> AsianCall(T=1,d=12,K=100,r=0.05,sigma=0.2,S0=100,method="naive",
+           sampling="MC",sampar=list(n=10^5))
          result error estimate
price 6.16783497   0.0527644567
delta 0.59392926   0.0032910230
gamma 0.03071819   0.0007318875
> AsianCall_AppLord(T = 1, d = 12, K = 100, r = 0.05, sigma = 0.25, S0 = 100, all = TRUE)
[1] 7.312461
> BS_EC(K=100, r = 0.05, sigma = 0.2, T = 0.25, S0 = 100)
     price      delta      gamma 
4.61499713 0.56946018 0.05694602 
> BS_EP(K=100, r = 0.05, sigma = 0.2, T = 0.25, S0 = 100)
      price       delta       gamma 
 3.37277718 -0.43053982  0.05694602 
> price = 1200 # current price of the bond
> C = 40 # coupon payment
> T= 30 # time to maturity
> par = 1000 # par value of the bond
> r = seq(0.02, 0.05, length = 300)