Hi, can anyone help me with the following code? Thanks!
library(mvtnorm)
f2 <- function(n, rho) {
var <- matrix(c(1,rho,rho,1), nrow=2, ncol=2, byrow=T)
beta <- seq(0, 1, length.out=n+1)
alpha <- sort (sapply(1-beta, qnorm))
x <- array(0, dim=c(n, n))
for (s in 1:n) {
for (t in 1:n){
if (s>=t)
x[s,t] <- pmvnorm(lower=c(alpha[s], alpha[t]),upper=c(alpha[s+1],
alpha[t+1]), mean=c(0,0), corr=var)/(s*(s+1))
else
x[s,t] <- pmvnorm(lower=c(alpha[s], alpha[t]),upper=c(alpha[s+1],
alpha[t+1]), mean-c(0,0), corr=var)/(t*(t+1))
} }
n*beta[2]-(n-1)*pmvnorm(lower=c(alpha[n-1],alpha[n-1]), upper=rep(Inf,2),
corr=var, mean=c(0,0))+n*(n-1)*sum(x) }
I got the error informatio as below if I want to calculate the value the
function for some specific parameters.
> f2(5, 0.1)
Error in checkmvArgs(lower = lower, upper = upper, mean = mean, corr corr, :
‘diag(corr)’ and ‘lower’ are of different length
[[alternative HTML version deleted]]
Did you try using:debug(f2) and then running the function, and see what is your variables condition before the error ? On Sun, Apr 19, 2009 at 12:58 PM, li li <hannah.hlx@gmail.com> wrote:> Hi, can anyone help me with the following code? Thanks! > > library(mvtnorm) > f2 <- function(n, rho) { > var <- matrix(c(1,rho,rho,1), nrow=2, ncol=2, byrow=T) > beta <- seq(0, 1, length.out=n+1) > alpha <- sort (sapply(1-beta, qnorm)) > x <- array(0, dim=c(n, n)) > for (s in 1:n) { > for (t in 1:n){ > if (s>=t) > x[s,t] <- pmvnorm(lower=c(alpha[s], alpha[t]),upper=c(alpha[s+1], > alpha[t+1]), mean=c(0,0), corr=var)/(s*(s+1)) > else > x[s,t] <- pmvnorm(lower=c(alpha[s], alpha[t]),upper=c(alpha[s+1], > alpha[t+1]), mean-c(0,0), corr=var)/(t*(t+1)) > } } > > n*beta[2]-(n-1)*pmvnorm(lower=c(alpha[n-1],alpha[n-1]), upper=rep(Inf,2), > corr=var, mean=c(0,0))+n*(n-1)*sum(x) } > > I got the error informatio as below if I want to calculate the value the > function for some specific parameters. > > > f2(5, 0.1) > Error in checkmvArgs(lower = lower, upper = upper, mean = mean, corr > corr, : > ‘diag(corr)’ and ‘lower’ are of different length > > [[alternative HTML version deleted]] > > > ______________________________________________ > R-help@r-project.org mailing list > stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide > R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. > >-- ---------------------------------------------- My contact information: Tal Galili Phone number: 972-50-3373767 FaceBook: Tal Galili My Blogs: r-statistics.com talgalili.com biostatistics.co.il [[alternative HTML version deleted]]