Correcy way to calculate positive NPV after monte carlo simulation

I am wondering if anyone can provide me links or idea about how i can calculate stochastic npv after monte carlo simulation and how i can calculate probability of npv>0? We first calculated deterministic npv with all the assumptions and then I took some important parameters where I can assign uncertainty and then assigned them uniform distribution (runif). But the probability of positive npv seems to be 0/1, nowhere in between, Is there something wrong with how i am calculating probability of positive npv?or how i am calculating npv_vec[i]?

...

ndraw <- 10000
set.seed(12345)
a_vec<- runif(100,10,20)
b_vec<- runif(100,20,30)
npv_vec <- rep(NA,ndraw)
profit_vec <- rep(NA,ndraw)

for(i in 1:ndraw) {
npv_vec[i] <- NPV_fun(a_vec[i],b_vec[i])
profit_vec[i] <- ifelse(npv_vec[i]>0,1,0)
}

calculate the probability of positive npv

pb_profit <- mean(profit_vec)
pb_profit

...

Is pb_profit is the right way to calculate that npv>0?

a_vec and b_vec don't seem to be defined in the code you've posted.

Thank you for pointing this out. I corrected them. a_vec and b_vec are basically some parameters (lets say capital cost and operational cost) that I consider to have uncertainty. a and b are already modeled in deterministic analysis, but now for monte carlo simulation, i am defining them as a_vec and b_vec. I am just wondering whether we should calculate the pb_proift from profit_vec[I] or pb_proift<- mean( profit_vec) is right to calculate positive npv? The code "pb_profit" is giving me something (o or 1), but i am just not sure that this is actually right way to calculate positive npv after monte carlo simulation

Every element of profit_vec is going to be a zero or one. If pb_profit equals exactly 0 or 1, that means that all the npv are positive or all or negative.

What's the NPV_fun() function?

When I include many uncertain variables (I incorporate at least 15/20 variables) in the for loop function, pb_profit sometimes gives me .75/.90 even. So it is perhaps to do with how many variables I am incorporating in uncertainty. The main npv function is as follows:

...

NPV_fun <- function(p_jet=24,capex=23,
opex=27) #default value
{

cash flow for Y0: initial investment

CF0 <- capex0.08(1+r_dis)^2+capex0.6(1+r_dis)+capex*0.32

cash flow for Y1-Y20: sale-expenses

exp_ut <- q_waterp_water+q_elecp_elec+
q_gasp_gas+q_h2p_h2
exp <- exp_ut+q_fdp_fd+q_labp_lab #275.56

print(exp)

sale_by <- q_dslp_dslp_oil+q_lpgp_lpgp_oil+q_napp_napp_oil+
q_mealp_meal
sale <- sale_by+q_jet
p_jet*p_oil

cat("revenue share for meal:",q_meal*p_meal/sale,"\n")

ebitda <- sale-exp

dep <- capex*0.8/t_dep
ebit <- ebitda-dep

other cash flow for Y1-10: mortgage, (depreciation)

amort <- amort.table(Loan=capex*debt,n=t_loan,pmt=NA,r_dis)[["Schedule"]]
CF1 <- as.data.frame(amort)
CF1$Depreciation <- dep

other cash flow for Y11-20: no mortgage, no depreciation

CF2 <- CF1*0

CF <- rbind(CF1,CF2)
CF <- cbind(Year=1:nrow(CF),CF)
CF$EBIT <- ebit
CF$EBT <- CF$EBIT-CF$Interest Paid
CF$Flow <- CF$EBT*(1-r_tax)+CF$Depreciation-CF$Principal Paid
CF$Discount <- r_dis
CF$PV <- CF$Flow/(1+CF$Discount)^CF$Year
NPV <- sum(CF$PV)-CF0

print(NPV)

}

NPV_base <- NPV_fun()

#######################################

Monte Carlo simulation

#######################################

ndraw <- 10000
set.seed(12345)
p_jet_vec<- runif(10000,10,20)
capex_vec<- runif(10000,20,30)
opex_vec<- runif(10000,30,40)

npv_vec <- rep(NA,ndraw)
profit_vec <- rep(NA,ndraw)

for(i in 1:ndraw) {
npv_vec[i] <- NPV_fun(p_jet_vec[i],capex_vec[i]
opex_vec[i])
profit_vec[i] <- ifelse(npv_vec[i]>0,1,0)
}

pb_profit <- mean(profit_vec)
pb_profit

...

As you can see, I put some random no, because I cant put original no here, which is sensitive. Hope this is good enough to understand?

You're running 10,000 Monte Carlo draws, but the three variables you define as inputs are only 100 long.

Sorry, it is corrected now. As I said, no coding mistake or issues occur during the running. Its just that I am concerned whether this is the right way to calculate probability of positive npv. every debuging issues has been already taken care of, so if you see a mistake in naming the function or anything, its just the mistake here, in my R code, everything runs smoothly

Yep, that is the right way to calculate the probability of a positive NPV.

Is there any websites/R book where I can take guidance from regarding this calculation? because i am not sure why, if i introduce some variables, it gives me p values such as .7/.8 but if i introduce say 14 variables, it gives me 100% probability, but then i introduce 23 variables into the for loop function for npv_vec[i], it gives me 0% probability of gettting positive npv. i dont understand why so much fluctuations please?

Also, in literatures, they usually provide mean, standard deviation, minimum and maximum values of this stochastic npv after monte carlo simulation. Am i to calculate these values from npv_vec or npv_vec[i] please?

One easy thing you might try is using hist() to draw histograms of NPV in various cases to see if that reveals anything.

from npv_vec[i], i get only one value as min, max, mean, standard deviaton. but from npv_vec, i get different vaules for all of these.

npv_vec[i] is one number, the value depending on i. npv_vec is the collection of all 10,000 values of npv_vec[i].

so for mean, max and min, i should try describe(npv_vec), not npv_vec[i]?

Exactly, you got it!

I have only skimmed this thread, but I have 2 general comments:

When you include a lot of variables you are going to get many that test at a level of significance. In other words, if you include 100 variables that are complete nonsense and unrelated, you will get 5 that will test as significant at a 95% level of confidence, and 10 at the 90%. You are far better offf to pick some of the more important, independent variables and test those.

This appears to be a problem that could benefit from a Bayesian approach where your data is used to test the values of the different prior parameters using a something like a JAGS program. That way, as you get new data, you can get better and better results.

Also, the posterior parameters are in the form of distributions, not points, so you can get answers like what is the probablility that the NPV is greater than some X, for example.

While this may far more work than is justified for your project, I would highy recommend John Kruschke's book: Doing Bayesian Data Analysis.

But, again, this could only be justified if the decisions made from your work were really important.