non è stato possibile raggiungere la piena precisione in 'qbeta'

Question

Sto eseguendo la versione R 2.14.0 su un PC che utilizza Windows 7 Ultimate (processore Intel Core i5-2400 3GHz con RAM da 8,00 GB). Fammi sapere se sono necessarie altre specifiche.

Sto provando a simulare i dati distribuiti beta correlati. Il metodo che sto usando è un'estensione di ciò che è scritto in questo documento:

http://onlinelibrary.wiley.com/doi/10.1002/asmb.901/pdf

1. In sostanza, comincio simulando dati normale multivariata (usando mvrnorm() funzione dal MASS).
2. Quindi io uso pnorm() per applicare la trasformazione probit a questi dati in modo che i miei nuovi dati vettoriali vivano su (0,1). E sono ancora correlati secondo la dichiarazione precedente.
3. Quindi, dati questi dati trasformati probit, applico la funzione qbeta() con determinati parametri shape1 e shape2, per recuperare i dati beta correlati con determinate proprietà di media e dispersione.

Conosco altri metodi per generare dati beta correlati. Sono interessato al motivo per cui lo qbeta() causa il fallimento di questo metodo per alcuni "semi". Di seguito è il messaggio di errore che ottengo.

Warning message: 
In qbeta(probit_y0, shape1 = a0, shape2 = b0) : 
    full precision may not have been achieved in 'qbeta' 

Cosa significa? Come può essere evitato? Quando si verifica nel contesto di una simulazione più ampia, qual è il modo migliore per garantire che questo problema non risolva l'intero codice di simulazione di origine (utilizzando source())?

Ho eseguito il seguente codice per i semi interi da 1: 1000. Seed = 899 era l'unico valore che mi dava problemi. Anche se è problematico qui, inevitabilmente sarà problematico anche per altri semi.

library(MASS) 
set.seed(899) 
n0 <- 25 
n1 <- 25  
a0 <- 0.25  
b0 <- 4.75  
a1 <- 0.25  
b1 <- 4.75  
varcov_mat <- matrix(rep(0.25,n0*n0),ncol=n0) 
diag(varcov_mat) <- 1 
y0 <- mvrnorm(1,mu=rep(0,n0),Sigma=varcov_mat) 
y1 <- mvrnorm(1,mu=rep(0,n1),Sigma=varcov_mat) 
probit_y0 <- pnorm(y0) 
probit_y1 <- pnorm(y1) 
beta_y0 <- qbeta(probit_y0, shape1=a0, shape2=b0) 
beta_y1 <- qbeta(probit_y1, shape1=a1, shape2=b1) 

Il codice sopra riportato è un frammento di un pezzo di un progetto di simulazione più ampio. Ma il messaggio di avvertimento qbeta() è quello che mi sta facendo venire il mal di testa ora.

Qualsiasi aiuto che il gruppo potrebbe fornire sarebbe molto apprezzato.

Acclamazioni Chris

Answer 1

la ragione per l'errore è che l'algoritmo utilizzato per il calcolo Qbeta non convergono per quei valori dei parametri.

R utilizza AS 109 per calcolare qbeta (Cran, GW, KJ Martin e GE Thomas (1977). Nota AS R19 e Algorithm AS 109, Statistica applicata, 26, 111-114 e commenti successivi (AS83 e correzione) .). R tenta di calcolare il valore in 1000 iterazioni. Se non è possibile nelle 1000 iterazioni, ottieni il messaggio di errore che hai visto.

Questo è il codice per qbeta.