R問題中的基本蒙地卡羅現值計算
我正在自學蒙特卡羅應用程序以及對現值的應用程序。
但是,我使用的值是具有預定義最小值和最大值的均勻分佈變數。
我正在嘗試將中心極限定理應用於它,因此它應該接近正態分佈,從而為我提供比其他值更有可能的值。
現在我在將 CLT 應用到我目前的設置時遇到了一些困難。
#Empty variable to store the list in results = NULL #loop for(i in 1:1000){ #rate between 1% and 20% r <- runif(1, 0.01, 0.2) #expected capital return k <- runif(1, 500, 900) #periods p <- 2 #do the actual calculation pv <- k * (1 + r)^p #Store the calculation results <- rbind(results, pv) } histogram <- hist(results) plot(histogram) summary <- summary(results) print(summary) standarddeviation <- sd(results) print(standarddeviation) 165.108 V1 Min. : 512.2 1st Qu.: 727.8 Median : 848.8 Mean : 856.8 3rd Qu.: 974.8 Max. :1292.6
所以我從中得到的是平均值是 856.8,這就是該項目目前或多或少值得投資的項目。
我想知道我的方法論和推理是否正確,因為我覺得我可能已經走深了。
您所做的並不是真正應用中心極限定理 (CLT),而是應用大數定律。如果我正確理解了您的問題,您可以從以下資訊開始:
- 未來貼現率服從均勻分佈: $ r \sim U(1%, 20%) $ .
- 未來資本分佈均勻: $ k \sim U(500, 900) $ .
- $ r $ 和 $ k $ 是獨立的。
由此您可以輕鬆計算“真實”的預期現值: $$ \begin{align} \mathbb{E}[pv] &= \mathbb{E}[k \cdot (1+r)^p] = \mathbb{E}[k] \cdot \mathbb{E}[(1+r)^p] \[2mm] &= 700 \cdot \int_{0.01}^{0.2} (1+x)^p \frac {1}{0.2 - 0.01} ; dx \[2mm] &= 700 \cdot \frac{1.2^{p+1} - 1.01^{p+1}}{0.19 \cdot (p+1)}. \end{align} $$
為了 $ p = 2 $ 你明白了 $ \mathbb{E}[pv] = 856.8233 $ .
您的程式碼表明大數定律有效:如果您重複實驗的次數足夠多,那麼該實驗的平均結果將接近其“真實”預期值。
關於程式碼本身的最後一條評論。在使用隨機數時,建議在 R 中使用隨機種子。這可以使用
set.seed
函式來實現。通過設置隨機種子,您可以使您的結果具有可重複性。您也可以使用replicate
而不是 for 循環來使您的程式碼更具可讀性:my_pv <- function(p){ r <- runif(1, 0.01, 0.2) k <- runif(1, 500, 900) return(k * (1 + r)^p) } set.seed(1234) result <- replicate(1000, my_pv(2)) mean(result) [1] 854.8702