Why Monte Carlo code results in ONE despite no. of trials 100 million

[Abdo]

I'm getting ONE in Monte Carlo simulation, I increased the number of simulations till 1e8, but still ones I`m trying to obtain the probability of failure. Here is my code Inputs are the mean and standard deviation and the type of distribution of the three parameters shown in function z. Here is the Code in R-studio

#MCcrude<-function(M,S,P=NULL,DS=NULL,Nsm,Nrp){
Nvr=length(M)
pf=NULL
for(k1 in seq(1:Nrp)){
X=NULL
for(j1 in seq(1,Nvr)){
X=cbind(X,rnorm(Nsm,mean= M[j1],sd =S[j1]))
}
z=30*X[,2]*X[,1]-X[,2]X[,1]^2-30X[,3]*X[,1]-X[,3]*X[,1]
pf=rbind(pf,length(which(z<=0))/Nsm) #probability of failure
}
return(pf)
}#

[fredoxvii]

Hello,

I found an error in you z function. The code below is working for me.

M <- rnorm(n = 10, mean = 10, sd = 1)
S <- rnorm(n = 10, mean = 1, sd = 1)
Nrp <- 10
Nvr <- 10
Nsm <- 10

MCcrude<-function(M,S,P=NULL,DS=NULL,Nsm,Nrp){
Nvr=length(M)
pf=NULL
#pf <- numeric(Nrp)
for(k1 in seq(1:Nrp)){
  X=NULL
  #X <- numeric(Nrp)
  for(j1 in 1:Nvr){
    X=cbind(X,rnorm(Nsm,mean= M[j1],sd = S[j1]))
    #message(X)
    #X[j1] <- rnorm(Nsm,mean= M[j1],sd =S[j1])
  }
  z <- 30*X[,2]*X[,1]-X[,2]*X[,1]^2-30*X[,3]*X[,1]-X[,3]*X[,1] # error
  pf=rbind(pf,length(which(z<=0))/Nsm) #probability of failure
  #pf[Nrp] <- sum(z<=0)/Nsm
}
return(pf)
}

MCcrude(M = M, S = S, Nsm = 20, Nrp = 20)

Results:

> MCcrude(M = M, S = S, Nsm = 20, Nrp = 20)
      [,1]
 [1,] 0.85
 [2,] 0.95
 [3,] 0.95
 [4,] 0.95
 [5,] 1.00
 [6,] 1.00
 [7,] 0.90
 [8,] 0.95
 [9,] 0.95
[10,] 1.00
[11,] 0.95
[12,] 0.95
[13,] 0.95
[14,] 0.90
[15,] 0.95
[16,] 0.95
[17,] 0.90
[18,] 0.90
[19,] 1.00
[20,] 1.00
There were 50 or more warnings (use warnings() to see the first 50)

[fredoxvii]

I would also suggest that you not use Nvr more than 3 because you are creating more columns that the z function is referencing, thus, it is just wasted calculation time. With an Nvr of 10 I had 10 columns, but the max X referenced is X[,3]. Just a thought.

[Abdo]

@fredoxvii

M <- rnorm(n = 10, mean = 10, sd = 1)
S <- rnorm(n = 10, mean = 1, sd = 1)

Why you are assigning S and M as long vectors like this

```S <- rnorm(n = 10, mean = 1, sd = 1)
I have a vector of 3 elements and I assign the like this
```S<-c(15,20,4)
```M<-c(14,400,1000)

Thanks

[fredoxvii]

Hello,

Your original question didn't clarify the inputs, so I was testing different scenarios to be thorough.