Dear All
Below are calls to functions to calculate bivariate and univariate
logistic probabilities.It works for the following sample program (with
results p1=p2 and p3=p4), but similar calls in a more elaborated program
produced unpredicted results.
My question is whether I am doing something bad (which I should avoid)
in my calls to mycdf2 and mycdf to obtain p2 and p3, respectively. Thank
you.
Steven Yen
pbivlogis <- function(x,y,rho){
# *********************************************
# Bivariate logistic CDF
# *********************************************
p<-(1+exp(-x)+exp(-y)+(1-rho)*exp(-x-y))^(-1)
return(p)
}
mycdf <- function(q,logistic=FALSE){
# *********************************************
# Univariate CDF: normal or logistic
# *********************************************
if(!logistic){
p<-pnorm(q)
} else {
p<-plogis(q)
}
return(p)
}
mycdf2 <- function(x,y,rho,logistic=FALSE){
# *********************************************
# Calling bivariate CDF: normal or logistic
# *********************************************
if(!logistic){
p<-pbivnorm(x,y,rho,recycle=T)
} else {
p<-pbivlogis(x,y,rho)
}
return(p)
}
set.seed(123)
x<-runif(n=5,min=-3,max=3)
y<-runif(n=5,min=-2,max=4)
rho<-0.5
p1<-pbivlogis(x,y,rho); p1
p2<-mycdf2(x,y,rho,logistic=TRUE); p2
p3<-mycdf(x,logistic=T); p3
p4<-plogis(x); p4
Results
> set.seed(123)
> x<-runif(n=5,min=-3,max=3)
> y<-runif(n=5,min=-2,max=4)
> rho<-0.5
> p1<-pbivlogis(x,y,rho); p1
[1] 0.04937376 0.65977865 0.35821101 0.72243120 0.63881214
> p2<-mycdf2(x,y,rho,logistic=TRUE); p2
[1] 0.04937376 0.65977865 0.35821101 0.72243120 0.63881214
> p3<-mycdf(x,logistic=T); p3
[1] 0.2184819 0.8493908 0.3667608 0.9087199 0.9335661
> p4<-plogis(x); p4
[1] 0.2184819 0.8493908 0.3667608 0.9087199 0.9335661
>
______________________________________________
R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
