[Gsdesign-commits] r254 - in pkg/gsDesign: . R man

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Mon Jan 11 04:43:30 CET 2010


Author: keaven
Date: 2010-01-11 04:43:30 +0100 (Mon, 11 Jan 2010)
New Revision: 254

Modified:
   pkg/gsDesign/DESCRIPTION
   pkg/gsDesign/R/gsCP.R
   pkg/gsDesign/R/gsDesign.R
   pkg/gsDesign/R/gsqplot.R
   pkg/gsDesign/man/gsCP.Rd
   pkg/gsDesign/man/gsDesign.Rd
Log:
2.2-7 updates to R and man

Modified: pkg/gsDesign/DESCRIPTION
===================================================================
--- pkg/gsDesign/DESCRIPTION	2010-01-10 23:50:43 UTC (rev 253)
+++ pkg/gsDesign/DESCRIPTION	2010-01-11 03:43:30 UTC (rev 254)
@@ -1,5 +1,5 @@
 Package: gsDesign
-Version: 2.2-5
+Version: 2.2-7
 Title: Group Sequential Design
 Author: Keaven Anderson 
 Maintainer: Keaven Anderson <keaven_anderson at merck.com>

Modified: pkg/gsDesign/R/gsCP.R
===================================================================
--- pkg/gsDesign/R/gsCP.R	2010-01-10 23:50:43 UTC (rev 253)
+++ pkg/gsDesign/R/gsCP.R	2010-01-11 03:43:30 UTC (rev 254)
@@ -50,9 +50,9 @@
         stop("gsCP must have x$lower$bound[i]<=zi<=x$upper$bound[i]")            
     }
   
-    if (!is.real(theta) || is.na(theta))
+    if (is.null(theta))
     {
-        theta <- c(zi/sqrt(x$n.I[i]), x$theta)
+        theta <- c(zi/sqrt(x$n.I[i]), 0, x$delta)
     }
     
     knew <- x$k-i
@@ -76,56 +76,19 @@
     {    
         stop("gsPP must be called with class of x either gsProbability or gsDesign")
     }
-    
-    if (i < 1 || i >= x$k)
-    {    
-        stop("gsPP must be called with i from 1 to x$k-1")
-    }
-    
     test.type <- ifelse(is(x, "gsProbability"), 3, x$test.type)
-    
-    if (zi > x$upper$bound[i])
-    {    
-        stop("gsPP must have x$lower$bound[i] <= zi <= x$upper$bound[i]")
-    }
-    else if (test.type > 1 && zi < x$lower$bound[i])
-    {
-        stop("gsPP must have x$lower$bound[i]<=zi<=x$upper$bound[i]")            
-    }
-  
-    if (!is.real(theta) || is.na(theta))
-    {
-        theta <- c(zi/sqrt(x$n.I[i]), x$theta)
-    }
-    
-    if (i < 1 || i >= x$k)
-    {    
-        stop("gsPP must be called with i from 1 to x$k-1")
-    }
-    
-    test.type <- ifelse(is(x, "gsProbability"), 3, x$test.type)
-    
-    if (zi > x$upper$bound[i])
-    {    
-        stop("gsCP must have x$lower$bound[i] <= zi <= x$upper$bound[i]")
-    }
-    else if (test.type > 1 && zi < x$lower$bound[i])
-    {
-        stop("gsCP must have x$lower$bound[i]<=zi<=x$upper$bound[i]")            
-    }
-  
-    if (!is.real(theta) || is.na(theta))
-    {
-        theta <- c(zi/sqrt(x$n.I[i]), x$theta)
-    }
+    checkScalar(i, "integer", c(1, x$k-1))
+    checkScalar(zi, "numeric", c(ifelse(test.type==1, -20, x$lower$bound[i]), x$upper$bound[i]), c(TRUE,TRUE))
+    checkVector(wgts, "numeric", c(-Inf, Inf), c(FALSE, FALSE))
+    checkVector(theta, "numeric", c(-Inf,Inf), c(FALSE, FALSE))
+    checkLengths(theta, wgts)
+    checkScalar(r, "integer", c(1, 80))
 
     cp <- gsCP(x = x, i = i, theta = theta, zi = zi, r = r)
-    cp$upper$prob %*% wgts   
     gsDen <- gsDensity(x, theta=theta, i=i, zi=zi, r=r)
-    if (!total) return(cp$upper$prob %*% t(gsDen$density))
-    one <- array(1, x$k-i)
-    cpmarg <- one %*% cp$upper$prob
-    sum(gsDen$density * cpmarg)
+	 pp <- cp$upper$prob %*% t(gsDen$density * wgts) / sum(gsDen$density * wgts)
+    if (total) return(sum(pp))
+    else return(pp)
 }
 
 "gsBoundCP" <- function(x, theta="thetahat", r=18)

Modified: pkg/gsDesign/R/gsDesign.R
===================================================================
--- pkg/gsDesign/R/gsDesign.R	2010-01-10 23:50:43 UTC (rev 253)
+++ pkg/gsDesign/R/gsDesign.R	2010-01-11 03:43:30 UTC (rev 254)
@@ -120,14 +120,14 @@
 
 "gsDesign"<-function(k=3, test.type=4, alpha=0.025, beta=0.1, astar=0,  
         delta=0, n.fix=1, timing=1, sfu=sfHSD, sfupar=-4,
-        sfl=sfHSD, sflpar=-2, tol=0.000001, r=18, n.I=0, maxn.IPlan=0, nFixSurv=0) 
+        sfl=sfHSD, sflpar=-2, tol=0.000001, r=18, n.I=0, maxn.IPlan=0, nFixSurv=0, endpoint=NULL) 
 {
     # Derive a group sequential design and return in a gsDesign structure
     
     # set up class variable x for gsDesign being requested
     x <- list(k=k, test.type=test.type, alpha=alpha, beta=beta, astar=astar,
             delta=delta, n.fix=n.fix, timing=timing, tol=tol, r=r, n.I=n.I, maxn.IPlan=maxn.IPlan,
-            nFixSurv=nFixSurv, nSurv=0)
+            nFixSurv=nFixSurv, nSurv=0, endpoint=endpoint)
     
     class(x) <- "gsDesign"
     
@@ -991,7 +991,7 @@
         {
             stop("maxn.IPlan can only be > 0 if spending functions are used for boundaries")
         }
-        
+
         x$timing <- x$n.I / x$maxn.IPlan
         
         if (x$n.I[x$k-1] >= x$maxn.IPlan)

Modified: pkg/gsDesign/R/gsqplot.R
===================================================================
--- pkg/gsDesign/R/gsqplot.R	2010-01-10 23:50:43 UTC (rev 253)
+++ pkg/gsDesign/R/gsqplot.R	2010-01-11 03:43:30 UTC (rev 254)
@@ -174,8 +174,8 @@
 
 "plotgsCP" <- function(x, theta="thetahat", main="Conditional power at interim stopping boundaries", 
         ylab=NULL, geom="line",
-        xlab=ifelse(x$n.I[x$k] < 3, "Sample size relative to fixed design", "N"), xlim=NULL,
-        lty=1, col=1, pch=22, textcex=1, legtext=gsLegendText(test.type), dgt=3, nlabel=TRUE, 
+        xlab=ifelse(x$n.fix == 1, "Sample size relative to fixed design", "N"), xlim=NULL,
+        lty=1, col=1, pch=22, textcex=1, legtext=NULL,  dgt=3, nlabel=TRUE, 
         base=FALSE, ...)
 {    
 	if (is.null(ylab))
@@ -187,8 +187,9 @@
 	if (!is.numeric(xlim))
 	{	xlim <- range(x$n.I[1:(x$k-1)])
 		xlim <- xlim + c(-.05,.05) * (xlim[2] - xlim[1])
+		if (x$k==2) xlim=xlim+c(-1, 1)
 	}
-	if(x$n.I[x$k] < 3) 
+	if(x$n.fix == 1) 
 	{	nround <- 3 
 		ntx <- "r="
 		if (is.null(xlab)) xlab <- "Information relative to fixed sample design"
@@ -198,6 +199,7 @@
 		if (is.null(xlab)) xlab <- "N"
 	}
 	test.type <- ifelse(is(x,"gsProbability"), 3, x$test.type)    
+	if (is.null(legtext)) legtext <- gsLegendText(test.type)
 	y <- gsBoundCP(x, theta=theta)
 	ymax <- 1.05
 	ymin <- - 0.1
@@ -215,7 +217,8 @@
     
     if (test.type > 1)
     {
-        ymid <- (y[2, 2] + y[2, 1]) / 2
+        if (x$k > 2) ymid <- (y[2, 2] + y[2, 1]) / 2
+        else ymid <- mean(y)
     }
     
 	if (base)
@@ -235,7 +238,6 @@
 		}
 		text(xtext, ymid, legtext[2], cex=textcex)
 		text(xtext, 1.03, legtext[1], cex=textcex)
-		invisible(x)
 	}else
 	{	N <- as.numeric(x$n.I[1:(x$k-1)])
 		CP <- y[,2] 
@@ -261,7 +263,7 @@
 		if (base)
 		{	#text(x=x$n.I[1:(x$k-1)], y=array(ymin/2, x$k-1), as.character(round(x$n.I[1:(x$k-1)],nround)), cex=textcex)
 		}
-		if (max(x$n.I) < 3)
+		if (x$n.fix == 1)
 		{	if (base)
 			{	text(x=y2$N, y=y2$CP, paste(array("r=",x$k), y2$Ztxt, sep=""), cex=textcex)
 			}else
@@ -269,7 +271,7 @@
 			}
 		}else
 		{	if(base)
-			{	text(x=y2$N, y=y2$Z, paste(array("N=",x$k), y2$Ztxt, sep=""), cex=textcex)
+			{	text(x=y2$N, y=y2$CP, paste(array("N=",x$k), y2$Ztxt, sep=""), cex=textcex)
 			}else
 			{	p <- p + geom_text(data=y2,label=paste(array("N=",x$k-1), y2$Ztxt, sep=""))
 			}

Modified: pkg/gsDesign/man/gsCP.Rd
===================================================================
--- pkg/gsDesign/man/gsCP.Rd	2010-01-10 23:50:43 UTC (rev 253)
+++ pkg/gsDesign/man/gsCP.Rd	2010-01-11 03:43:30 UTC (rev 254)
@@ -5,9 +5,8 @@
 \alias{gsCPOS}
 \title{2.4: Conditional and Predictive Power, Overall and Conditional Probability of Success}
 \description{
-\code{gsCP()} takes a given group sequential design, assumes an interim z-statistic 
-at a specified interim analysis and computes boundary crossing probabilities at future planned analyses.
-\code{gsPP()} averages conditional power across a prior distribution to produce predictive power.
+\code{gsCP()} computes conditional boundary crossing probabilities at future planned analyses for a given group sequential design assuming an interim z-statistic at a specified interim analysis.
+\code{gsPP()} averages conditional power across a prior distribution to compute predictive power.
 \code{gsPOS()} computes the probability of success for a trial using a prior distribution to average power over a set of \code{theta} values of interest.
 \code{gsCPOS()} assumes no boundary has been crossed before and including an interim analysis of interest, and computes the probability of success based on this event.
 Note that \code{gsCP()} and \code{gsPP()} take only the interim test statistic into account in computing conditional probabilities, while \code{gsCPOS()} conditions on not crossing any bound through a specified interim analysis.
@@ -33,15 +32,15 @@
 }
 \value{
 \code{gsCP()} returns an object of the class \code{gsProbability}. 
-Based on the input design and the interim test statistic, the output object has bounds for test statistics
-computed based on observations after interim \code{i} that are equivalent to the original design crossing boundaries conditional on the interim test statistic value input. 
+Based on the input design and the interim test statistic, the output gsDesign object has bounds for test statistics computed based on solely on observations after interim \code{i}. 
 Boundary crossing probabilities are computed for the input \eqn{\theta}{theta} values.
+See manual and examples.
 
-\code{gsPP()} returns a real value indicating the predictive power of the trial conditional on the interim test statistic \code{zi} at analysis \code{i}.
+\code{gsPP()} if total==TRUE, returns a real value indicating the predictive power of the trial conditional on the interim test statistic \code{zi} at analysis \code{i}; otherwise returns vector with predictive power for each future planned analysis.
 
 \code{gsPOS()} returns a real value indicating the probability of a positive study weighted by the prior distribution input for \code{theta}.
 
-\code{gsCPOS()} returns a real value indicating the probability of a positive study weighted by the prior distribution input for \code{theta} conditional on an interim test statistic \code{zi} at analysis \code{i}.
+\code{gsCPOS()} returns a real value indicating the probability of a positive study weighted by the prior distribution input for \code{theta} conditional on an interim test statistic \code{zi} in the continuation region at analysis \code{i}.
 }
 
 
@@ -64,17 +63,50 @@
 x <- gsDesign(k=5)
 x
 
-# assuming a z-value of .5 at analysis 2, what are conditional 
+# set up a prior distribution for the treatment effect
+# that is normal with mean .75*x$delta and standard deviation x$delta/2
+prior <- normalGrid(mu=.75*x$delta, sigma=x$delta/2)
+names(prior)
+
+# compute POS for the design given the above prior distribution for theta
+gsPOS(x=x, theta=prior$z, wgts=prior$wgts)
+
+# assume POS should only count cases in prior where theta >= x$delta/2
+gsPOS(x=x, theta=prior$z, wgts=prior$wgts*(prior$z>=x$delta/2))
+
+# assuming a z-value at lower bound at analysis 2, what are conditional 
 # boundary crossing probabilities for future analyses
 # assuming theta values from x as well as a value based on the interim
 # observed z
-CP <- gsCP(x, i=2, zi=.5)
+CP <- gsCP(x, i=2, zi=x$lower$bound[2])
 CP
 
 # summing values for crossing future upper bounds gives overall
 # conditional power for each theta value
 CP$theta
-CP$upper$prob %*% c(1,1,1)
+t(CP$upper$prob) %*% c(1,1,1)
+
+# compute predictive probability based on above assumptions
+gsPP(x, i=2, zi=x$lower$bound[2], theta=prior$z, wgts=prior$wgts)
+
+# if it is known that boundary not crossed at interim 2, use
+# gsCPOS to compute conditional POS based on this
+gsCPOS(x=x, i=2, theta=prior$z, wgts=prior$wgts)
 }
+
+# 2-stage example to compare results to direct computation
+x<-gsDesign(k=2)
+z1<- 0.5
+n1<-x$n.I[1]
+n2<-x$n.I[2]-x$n.I[1]
+thetahat<-z1/sqrt(n1)
+theta<-c(thetahat, 0 , x$delta)
+
+# conditional power direct computation
+pnorm((n2*theta+z1*sqrt(n1)-z*sqrt(n1+n2))/sqrt(n2))
+gsCP(x=x, zi=z1, i=1)$upper$prob
+
+
+
 \seealso{\code{\link{gsDesign}}, \code{\link{gsProbability}}, \code{\link{gsBoundCP}}}
 \keyword{design}

Modified: pkg/gsDesign/man/gsDesign.Rd
===================================================================
--- pkg/gsDesign/man/gsDesign.Rd	2010-01-10 23:50:43 UTC (rev 253)
+++ pkg/gsDesign/man/gsDesign.Rd	2010-01-11 03:43:30 UTC (rev 254)
@@ -7,7 +7,7 @@
 gsDesign(k=3, test.type=4, alpha=0.025, beta=0.1, astar=0,  
          delta=0, n.fix=1, timing=1, sfu=sfHSD, sfupar=-4,
          sfl=sfHSD, sflpar=-2, tol=0.000001, r=18, n.I = 0,
-         maxn.IPlan = 0, nFixSurv=0) 
+         maxn.IPlan = 0, nFixSurv=0, endpoint=NULL) 
 
 print.gsDesign(x,...)}
 
@@ -60,7 +60,8 @@
 	compute the total sample size required as well as the number of events at each analysis that will be returned in \code{n.fix}; this is rounded up to an even number.}|
 	\item{x}{In \code{print.gsDesign} this is an object of class gsDesign.}
 	\item{...}{This should allow optional arguments that are standard when calling \code{print}.}
-}
+\code{endpoint}{An optional character string that should represent the type of endpoint used for the study. This may be used by output functions. Types most likely to be recognized initially are "TTE" for time-to-event outcomes with fixed design sample size generated by \code{nSurvival()} and "Binomial" for 2-sample binomial outcomes with fixed design sample size generated by \code{nBinomial()}.
+}}
 \value{
 An object of the class \code{gsDesign}. This class has the following elements and upon return from 
 \code{gsDesign()} contains:
@@ -89,6 +90,7 @@
 If \code{delta=0} was input, then this is the sample size required for the specified group sequential design when a fixed design requires a sample size of \code{n.fix}.
 If \code{delta=0} and \code{n.fix=1} then this is the relative sample size compared to a fixed design; see details and examples.}
 \item{maxn.IPlan}{As input.}
+\item{endpoint}{As input.}
 }
 \details{
 Many parameters normally take on default values and thus do not require explicit specification.



More information about the Gsdesign-commits mailing list