Monte Carlo con R

Buongiorno, devo scrivere un algoritmo con il metodo di Monte Carlo, tale che riesca a calcolare il volume di una sfera di raggio=1 centrata in (0,0,0). Questo è l'algoritmo base:
# intMonteCarlo: funzione per il calcolo di integrali con il metodo Monte Carlo grezzo. L'integrale viene calcolato come media dei valori della funzione in punti campionati da distribuzione uniforme sull'intervallo di integrazione. La stima dell'errore e' una misura di quanto la funzione si discosta dal valore medio sull'intervallo di integrazione.
Input:
- fcn: funzione integranda
- a: limite inferiore integrazione
- b: limite superiore integrazione
- n: numero estrazioni
- analisi: se TRUE si ritorna stima integrale, e sigma per analisi con n crescente,
default FALSE
Output:
- media e sigma se analisi = FALSE
- dataframe con vettori N, integrale(N), sigma, se analisi = TRUE

intMonteCarlo <- function(fcn, a, b, n, iop = FALSE)
{
  x <- runif(n, min = a, max = b)
  fX <- sapply(x, fcn)
  if (iop)
  {
    N <- seq(1,n)
    stimaInt <- (b-a) * (cumsum(fX)/(1:n))
    stimaErr <- (b-a) * sqrt((cumsum(fX^2)/N-(cumsum(fX)/N)^2)/N)
    return(data.frame(N, stimaInt, stimaErr))
  } else
  {
    stimaInt <- (b-a) * mean(fX)
    stimaErr <- (b-a) * sqrt((sum(fX^2)/n - mean(fX)^2)/n)
    return(c(stimaInt, stimaErr))
  }
}

come faccio ad implementarlo per una funzione in 3 dimensioni, quale è quella della sfera(x^2 + y^2 + z^2 -1)?