R version 2.12.1 (2010-12-16) Copyright (C) 2010 The R Foundation for Statistical Computing ISBN 3-900051-07-0 Platform: i486-pc-linux-gnu (32-bit) R is free software and comes with ABSOLUTELY NO WARRANTY. You are welcome to redistribute it under certain conditions. Type 'license()' or 'licence()' for distribution details. R is a collaborative project with many contributors. Type 'contributors()' for more information and 'citation()' on how to cite R or R packages in publications. Type 'demo()' for some demos, 'help()' for on-line help, or 'help.start()' for an HTML browser interface to help. Type 'q()' to quit R. > x <- c(0.224565646 + ,0.261571480 + ,0.253029834 + ,0.033938060 + ,0.353208554 + ,0.033295045 + ,0.002171707 + ,0.084136132 + ,0.077468368 + ,0.010160325 + ,0.079746191 + ,0.263485597 + ,0.460389066 + ,0.132348066 + ,0.529154319 + ,0.032287005 + ,0.472190992 + ,0.037985170 + ,0.019795866 + ,0.053763887 + ,0.086771175 + ,0.143764004 + ,0.202780264 + ,0.006932185 + ,0.124637036 + ,0.116166293 + ,0.037363442 + ,0.416532177 + ,0.612491026 + ,0.012781084 + ,0.163992436 + ,0.212717925 + ,0.092409609 + ,0.150134816 + ,0.047430792 + ,0.022017328 + ,0.350262865 + ,0.013821921 + ,0.578796918 + ,0.030264300 + ,0.090627282 + ,0.029354226 + ,0.182663201 + ,0.372085552 + ,0.235142217 + ,0.311617946 + ,0.030548270 + ,0.068424450 + ,0.043757662 + ,0.238109173 + ,0.465951209 + ,0.127494709 + ,0.010810575 + ,0.230141979 + ,0.426687450 + ,0.175605499 + ,0.005595384 + ,0.203565169 + ,0.236888269 + ,0.153586896 + ,0.027445414 + ,0.026240968 + ,0.299075759 + ,0.171469487 + ,0.016355307 + ,0.220064442 + ,0.137054994 + ,0.034683689 + ,0.085767702 + ,0.119746590 + ,0.223347602 + ,0.026601623 + ,0.063485064 + ,0.467803972 + ,0.012650088 + ,0.048496183 + ,0.047170846 + ,0.043196485 + ,0.095073370 + ,0.538513979 + ,0.045034181 + ,0.427235474 + ,0.289835816 + ,0.150884924 + ,0.176887675 + ,0.369202741 + ,0.178461483 + ,0.077501113 + ,0.281105892 + ,0.051375864 + ,0.357311398 + ,0.561861189 + ,0.028342902 + ,0.326209843 + ,0.303994601 + ,0.020530542 + ,0.278139180 + ,0.582716856 + ,0.071720238 + ,0.237546024) > par4 = '6' > par3 = '2' > par2 = '3' > par1 = '0.1' > #'GNU S' R Code compiled by R2WASP v. 1.0.44 () > #Author: Prof. Dr. P. Wessa > #To cite this work: AUTHOR(S), (YEAR), YOUR SOFTWARE TITLE (vNUMBER) in Free Statistics Software (v$_version), Office for Research Development and Education, URL http://www.wessa.net/rwasp_YOURPAGE.wasp/ > #Source of accompanying publication: Office for Research, Development, and Education > #Technical description: > library(MASS) > PPCCBeta <- function(shape1, shape2, x) + { + x <- sort(x) + pp <- ppoints(x) + cor(qbeta(pp, shape1=shape1, shape2=shape2), x) + } > par1 <- as.numeric(par1) > par2 <- as.numeric(par2) > par3 <- as.numeric(par3) > par4 <- as.numeric(par4) > if (par1 < 0.1) par1 <- 0.1 > if (par1 > 10) par1 <- 10 > if (par2 < 0.1) par2 <- 0.1 > if (par2 > 10) par2 <- 10 > if (par3 < 0.1) par3 <- 0.1 > if (par3 > 10) par3 <- 10 > if (par4 < 0.1) par4 <- 0.1 > if (par4 > 10) par4 <- 10 > par1h <- par1*10 > par2h <- par2*10 > par3h <- par3*10 > par4h <- par4*10 > sortx <- sort(x) > c <- array(NA,dim=c(par2h,par4h)) > for (i in par1h:par2h) + { + for (j in par3h:par4h) + { + c[i,j] <- cor(qbeta(ppoints(x), shape1=i/10,shape2=j/10),sortx) + } + } > postscript(file="/var/www/rcomp/tmp/1xifb1307089643.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > filled.contour((par1h:par2h)/10,(par3h:par4h)/10,c[par1h:par2h,par3h:par4h],xlab='shape1',ylab='shape2',main='PPCC Contour Plot - Beta') > dev.off() null device 1 > xbar <- mean(x) > xvar <- var(x) > (a <- (xbar*(1-xbar)/xvar - 1)*xbar) [1] 0.823442 > (b <- (1-xbar)*a/xbar) [1] 3.752365 > (f<-fitdistr(x, 'beta',list(shape1=a,shape2=b))) shape1 shape2 0.8831678 4.0315009 (0.1091238) (0.6117423) > xlab <- paste('Beta(shape1=',round(f$estimate[[1]],2)) > xlab <- paste(xlab,', shape2=') > xlab <- paste(xlab,round(f$estimate[[2]],2)) > xlab <- paste(xlab,')') > postscript(file="/var/www/rcomp/tmp/2q8q41307089643.ps",horizontal=F,onefile=F,pagecentre=F,paper="special",width=8.3333333333333,height=5.5555555555556) > myser <- qbeta(ppoints(x), shape1=f$estimate[[1]], shape2=f$estimate[[2]]) > qqplot(myser, x, main='QQ plot (Beta)', xlab=xlab ) > grid() > dev.off() null device 1 > > #Note: the /var/www/rcomp/createtable file can be downloaded at http://www.wessa.net/cretab > load(file="/var/www/rcomp/createtable") > > a<-table.start() > a<-table.row.start(a) > a<-table.element(a,'Parameter',1,TRUE) > a<-table.element(a,'Estimated Value',1,TRUE) > a<-table.element(a,'Standard Deviation',1,TRUE) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,'shape1',header=TRUE) > a<-table.element(a,f$estimate[1]) > a<-table.element(a,f$sd[1]) > a<-table.row.end(a) > a<-table.row.start(a) > a<-table.element(a,'shape2',header=TRUE) > a<-table.element(a,f$estimate[2]) > a<-table.element(a,f$sd[2]) > a<-table.row.end(a) > a<-table.end(a) > table.save(a,file="/var/www/rcomp/tmp/3r4ag1307089643.tab") > > try(system("convert tmp/1xifb1307089643.ps tmp/1xifb1307089643.png",intern=TRUE)) character(0) > try(system("convert tmp/2q8q41307089643.ps tmp/2q8q41307089643.png",intern=TRUE)) character(0) > > > proc.time() user system elapsed 3.056 0.156 3.191