Nothing
#################################################################
#
# stsubpop.R
#
#######################
# stepp subpopulation #
#######################
# supporting routine to generate subpopulation based just on the coviarate and window
generate.all <- function(sp, win, covariate, coltype = NULL, coltrt = NULL, trts = NULL, minsubpops = NULL) {
r1 <- win@r1
r2 <- win@r2
e1 <- win@e1
e2 <- win@e2
if (win@type == "tail-oriented") {
zvals <- unique(sort(covariate))
if (min(r1) < min(zvals) | max(r2) > max(zvals)) stop("Cannot create tail-oriented window")
# nsubpop <- length(r1) + length(r2) + 1 ### [SV 27.02.2021: shouldn't it be length(r1)+length(r2)?]
nsubpop <- length(r1) +length(r2)
npatsub <- rep(0, nsubpop)
I0 <- rep(1, nsubpop)
I1 <- rep(length(zvals), nsubpop)
# if (length(r1) != 0) { ### [SV 27.02.2021: shouldn't it be > 1?]
if (length(r1) > 1) {
for (i in 1:length(r1)) {
npatsub[i] <- sum(covariate <= r1[i])
sel <- which(zvals <= r1[i])
I1[i] <- sel[length(sel)]
}
npatsub[length(r1) + 1] <- length(covariate)
}
# if (length(r2) != 0) { ### [SV 27.02.2021: shouldn't it be > 1?]
if (length(r2) > 1) {
npatsub[1] <- length(covariate)
for (i in 1:length(r2)) {
k <- i + length(r1)
npatsub[k] <- sum(covariate >= r2[i])
I0[k] <- which(zvals >= r2[i])[1]
}
}
} else if (win@type == "sliding") {
if (r2 >= length(covariate)) {
errmsg <- "patspop (r2) must be strictly smaller than the number of unique covariate values."
stop(errmsg)
}
zvals <- unique(sort(covariate))
absfreq <- as.numeric(table(covariate))
cumfreq <- cumsum(absfreq)
I0 <- rep(0, 1000) # max size of I0 and I1 is 1000; these limit the
I1 <- rep(0, 1000) # max no. of stepp subpopulation to 1000
I0[1] <- 1
I1[1] <- sum(cumfreq < r2) + 1
stopflag <- 0
nsubpop <- 2
while (stopflag == 0) {
indinf <- I0[nsubpop - 1] + 1
while ((cumfreq[I1[nsubpop - 1]] - cumfreq[indinf - 1]) > r1) {
indinf <- indinf + 1
}
I0[nsubpop] <- indinf
indsup <- I1[nsubpop - 1]
while (((cumfreq[indsup] - cumfreq[I0[nsubpop]] + absfreq[I0[nsubpop]]) < r2) && (stopflag == 0)) {
indsup <- indsup + 1
stopflag <- 1 * (indsup == length(zvals))
}
I1[nsubpop] <- indsup
nsubpop <- nsubpop + 1
}
nsubpop <- nsubpop - 1
npatsub <- rep(0, nsubpop)
npatsub[1] <- cumfreq[I1[1]]
for (i in 2:nsubpop) npatsub[i] <- cumfreq[I1[i]] - cumfreq[I0[i] - 1]
I0 <- I0[1:nsubpop]
I1 <- I1[1:nsubpop]
} else if (win@type == "sliding_events") {
#
# create 0/1 treatment assignment
#
txassign <- coltrt
for (i in 1:length(txassign)) {
if (txassign[i] != trts[1] & txassign[i] != trts[2]) {
txassign[i] <- NA
}
if (txassign[i] == trts[1]) {
txassign[i] <- 1
}
else if (txassign[i] == trts[2]) {
txassign[i] <- 0
}
}
#
# get the overlapping subgroups
#
# zvals is the unique values of the covariate that has been sorted at which there has been an event
zvals <- unique(sort(covariate))
# absfreq is the total number of events of that covariate value corresponding to zvals
absfreq <- as.numeric(table(covariate))
cumfreq <- cumsum(absfreq)
#
# Create vector of number of Type 1's recorded for each covariate
# value for each treatment type
#
type1valsTrt0 <- type1valsTrt1 <- rep(0, length(zvals))
for (i in 1:length(covariate)) {
if (coltype[i] == 1) {
covariateIndices <- which(zvals == covariate[i])
if (length(covariateIndices) != 1) {
stop("Non-unique covariate detected in zvals vector or value missing.")
}
covariateIndex <- covariateIndices[1]
#####NEW CODE
if (txassign[i] == 0) {
type1valsTrt0[covariateIndex] <- type1valsTrt0[covariateIndex] + 1
}
if (txassign[i] == 1) {
type1valsTrt1[covariateIndex] <- type1valsTrt1[covariateIndex] + 1
}
}
}
# Create cummulative graph of typevals
cumtypevalsTrt0 <- cumsum(type1valsTrt0)
cumtypevalsTrt1 <- cumsum(type1valsTrt1)
I0 <- I1 <- rep(0, length(zvals))
#
# sets the boundaries of the first subpopulation
#
I0[1] <- 1
#
# Determines the minimum upper boundary that satisfies the value of e2 for both treatment types
#
upperTrt0 <- sum(cumtypevalsTrt0 < e2) + 1
upperTrt1 <- sum(cumtypevalsTrt1 < e2) + 1
I1[1] <- max(upperTrt0, upperTrt1)
#
# Check to see if the upper boundary of the first subpop is greater than
# the size of the entire population. This indicates that there are
# insufficient events to satisfy E2.
#
if (I1[1] >= length(covariate)) {
# obsInitialSubpopCreationFailures <- obsInitialSubpopCreationFailures + 1
print(paste("Warning: Insufficient Type 1 Events found creating initial subpopulation (",
"Trt0=", cumtypevalsTrt0[length(cumtypevalsTrt0)],
" Trt1=", cumtypevalsTrt1[length(cumtypevalsTrt1)], ")", sep=""))
stop("Throwing out dataset.")
}
#
# Find the boundaries of each subpopulation
#
stopflag <- 0
nsubpop <- 2
while (stopflag == 0) {
#
# find the lower boundary of the next subpopulation
#
# indinf is the presumptive lower bound of the next subpopulation
indinf <- I0[nsubpop - 1] + 1
while (((cumtypevalsTrt0[I1[nsubpop - 1]] - cumtypevalsTrt0[indinf - 1]) > e1) ||
((cumtypevalsTrt1[I1[nsubpop - 1]] - cumtypevalsTrt1[indinf - 1]) > e1)) {
indinf <- indinf + 1
}
I0[nsubpop] <- indinf
indsup <- I1[nsubpop - 1]
#
# find the upper boundary of the next subpopulation
#
while ((((cumtypevalsTrt0[indsup] - cumtypevalsTrt0[I0[nsubpop]] + type1valsTrt0[I0[nsubpop]]) < e2) ||
((cumtypevalsTrt1[indsup] - cumtypevalsTrt1[I0[nsubpop]] + type1valsTrt1[I0[nsubpop]]) < e2)) &&
(stopflag == 0)) {
indsup <- indsup + 1
stopflag <- (indsup == length(zvals))
}
I1[nsubpop] <- indsup
# increment total number of subpopulations for next loop iteration
nsubpop <- nsubpop + 1
}
# decrement nsubpop to number of found subpopulations
nsubpop <- nsubpop - 1
#
# If the last subpopulation is too small for either treatment arm (< e2)
# include it in the previous subpopulation so that e2 is not violated
#
if (((cumtypevalsTrt0[I1[nsubpop]] - cumtypevalsTrt0[I0[nsubpop] - 1]) < e2) ||
((cumtypevalsTrt1[I1[nsubpop]] - cumtypevalsTrt1[I0[nsubpop] - 1]) < e2)) {
I1[nsubpop - 1] <- I1[nsubpop]
nsubpop <- nsubpop - 1
}
#
# Make sure there are enough subpopulations
#
if (nsubpop < minsubpops) {
# obsInsufficientSubpopsErrors <- obsInsufficientSubpopsErrors + 1
print(paste("Warning: too few subpopulations (", nsubpop, ")", sep = ""))
stop("Throwing out dataset.")
}
npatsub <- rep(0, nsubpop)
npatsub[1] <- cumfreq[I1[1]]
for (i in 2:nsubpop) npatsub[i] <- cumfreq[I1[i]] - cumfreq[I0[i] - 1]
neventsubTrt0 <- neventsubTrt1 <- rep(0, nsubpop)
neventsubTrt0[1] <- cumtypevalsTrt0[I1[1]]
neventsubTrt1[1] <- cumtypevalsTrt1[I1[1]]
for (i in 2:nsubpop) {
neventsubTrt0[i] <- cumtypevalsTrt0[I1[i]] - cumtypevalsTrt0[I0[i] - 1]
neventsubTrt1[i] <- cumtypevalsTrt1[I1[i]] - cumtypevalsTrt1[I0[i] - 1]
}
I0 <- I0[1:nsubpop]
I1 <- I1[1:nsubpop]
}
medians <- rep(0, nsubpop)
minc <- rep(NA, nsubpop)
maxc <- rep(NA, nsubpop)
npats <- length(covariate)
subpop <- matrix(rep(0, (npats*nsubpop)), ncol = nsubpop)
for (i in 1:nsubpop) {
subpop[, i] <- (covariate >= zvals[I0[i]]) * (covariate <= zvals[I1[i]])
medians[i] <- round((median(covariate[subpop[, i] == 1])), digits = 2)
minc[i] <- zvals[I0[i]]
maxc[i] <- zvals[I1[i]]
}
# update the object
sp@win <- win
sp@covar <- covariate
sp@nsubpop <- nsubpop
sp@subpop <- subpop
sp@npatsub <- npatsub
sp@medianz <- medians
sp@minc <- minc
sp@maxc <- maxc
if (sp@win@type == "sliding_events") {
sp@neventsubTrt0 <- neventsubTrt0
sp@neventsubTrt1 <- neventsubTrt1
} else {
sp@neventsubTrt0 <- NULL
sp@neventsubTrt1 <- NULL
}
sp@init <- TRUE
return(sp)
}
setClass("stsubpop",
representation(
win = "stwin", # stepp window object
covar = "numeric", # vector of covariate of interest (V)
nsubpop = "numeric", # number of subpopulation generated
subpop = "ANY", # matrix of subpopulations
npatsub = "numeric", # count of each subpopulation
medianz = "numeric", # median of V for each subpopulation
minc = "numeric", # minimum of V for each subpopulation, actual
maxc = "numeric", # maximum of V for each subpopulation, actual
neventsubTrt0 = "null_or_numeric",
neventsubTrt1 = "null_or_numeric",
init = "logical" # initialized
)
)
setMethod("initialize", "stsubpop",
function(.Object) {
.Object@init <- FALSE
return(.Object)
}
)
setValidity("stsubpop",
function(object) {
if (!is.numeric(object@covar) || length(object@covar) < object@win@r2) {
print("Invalid argument.")
return (FALSE)
}
return(TRUE)
}
)
setGeneric("generate", function(.Object, win, covariate, coltype, coltrt, trts, minsubpops)
standardGeneric("generate"))
setMethod("generate",
signature="stsubpop",
definition = function(.Object, win, covariate, coltype, coltrt, trts, minsubpops) {
if (missing(minsubpops)) {
minsubpops <- 2
}
.Object <- generate.all(.Object, win, covariate, coltype, coltrt, trts, minsubpops)
return(.Object)
}
)
setGeneric("merge", function(.Object, mergevector)
standardGeneric("merge"))
setMethod("merge",
signature = "stsubpop",
definition = function(.Object, mergevector) {
subp.new <- new("stsubpop")
if (.Object@nsubpop != length(mergevector)) stop("invalid merge vector length.")
if (.Object@nsubpop == 1) {
# if there is only one subpopulation left, there is nothing to merge
subpop.new <- .Object
} else {
v1 <- c(mergevector[1], mergevector[-length(mergevector)])
run <- c(1, which(v1 != mergevector), length(mergevector) + 1)
merge <- mergevector[1]
subpop.matrix <- NULL
nsubpop.new <- 0
for (i in 1:(length(run) - 1)) {
beg <- run[i]
end <- run[i + 1] - 1
if (merge & (beg != end)) {
subpop.matrix <- cbind(subpop.matrix, apply(.Object@subpop[, beg:end], 1, max))
nsubpop.new <- nsubpop.new + 1
} else {
subpop.matrix <- cbind(subpop.matrix, .Object@subpop[, beg:end])
nsubpop.new <- nsubpop.new + end - beg + 1
}
merge <- !merge
}
medianz <- rep(0, nsubpop.new)
minc <- rep(0, nsubpop.new)
maxc <- rep(0, nsubpop.new)
covariate <- .Object@covar
for (i in 1:nsubpop.new) {
subpop.cov <- covariate[subpop.matrix[, i] == 1]
medianz[i] <- round((median(subpop.cov)), digits = 2)
minc[i] <- min(subpop.cov)
maxc[i] <- max(subpop.cov)
}
subp.new@win <- .Object@win
subp.new@covar <- .Object@covar
subp.new@nsubpop <- nsubpop.new
subp.new@npatsub <- apply(subpop.matrix,2,sum)
subp.new@subpop <- subpop.matrix
subp.new@medianz <- medianz
subp.new@minc <- minc
subp.new@maxc <- maxc
subp.new@neventsubTrt0 <- .Object@neventsubTrt0
subp.new@neventsubTrt1 <- .Object@neventsubTrt1
subp.new@init <- TRUE
}
return(subp.new)
}
)
setGeneric("edge", function(.Object, j, side)
standardGeneric("edge"))
setMethod("edge",
signature = "stsubpop",
definition = function(.Object, j, side) {
if (side=="L" & j ==1) stop("invalid arg j.")
# sort the covariate value
covariate <- .Object@covar
Z <- sort(covariate)
ZU <- unique(Z)
if (side == "L") {
left.edge <- .Object@minc[j-1]
left.i <- (which(ZU == left.edge))[1]
rside <- which(ZU == .Object@maxc[j-1])
right.i <- rside[length(rside)]
n <- min(right.i - left.i + 1, length(ZU)-right.i+1)
start <- left.i
} else if (side == "R") {
right.edge <- .Object@maxc[j+1]
rside <- which(ZU == right.edge)
right.i <- rside[length(rside)]
left.i <- (which(ZU == .Object@minc[j+1]))[1]
n <- min(right.i-left.i+1,left.i)
start <- left.i-n+1
} else {
stop("unknown side.")
}
# create a special stepp subpopulation just for this edge
subpop.matrix <- matrix(rep(0, (length(covariate) * n)), ncol = n)
medianz <- rep(0, n)
minc <- rep(0, n)
maxc <- rep(0, n)
# print(paste(side, "j=",j, "start=",start, "n=",n,ZU[start], ZU[start+n-1]))
# print(ZU)
for (i in 1:n) {
# print(c(i, ZU[start + i - 1], ZU[start + i + n - 1]))
subpop.matrix[,i ] <- (covariate >= ZU[start + i - 1]) * (covariate <= ZU[start + i+n - 2])
subpop.cov <- covariate[subpop.matrix[, i] == 1]
medianz[i] <- round((median(subpop.cov)), digits = 2)
minc[i] <- min(subpop.cov)
maxc[i] <- max(subpop.cov)
}
subp.new <- new("stsubpop")
subp.new@win <- .Object@win
subp.new@covar <- .Object@covar
subp.new@nsubpop <- n
subp.new@subpop <- subpop.matrix
subp.new@npatsub <- apply(subpop.matrix,2,sum)
subp.new@medianz <- medianz
subp.new@minc <- minc
subp.new@maxc <- maxc
subp.new@neventsubTrt0 <- .Object@neventsubTrt0
subp.new@neventsubTrt1 <- .Object@neventsubTrt1
subp.new@init <- TRUE
return(subp.new)
}
)
setMethod("summary",
signature = "stsubpop",
definition = function(object) {
summary(object@win)
if (object@init) {
if (object@win@type == "sliding_events") {
write(paste("Number of subpopulations created:", object@nsubpop), file = "")
cat("\n")
write("Subpopulation summary information", file = "")
nper <- apply(object@subpop, 2, sum)
temp <- matrix(c(1:object@nsubpop, object@medianz, round(object@minc, digits = 4),
round(object@maxc, digits = 4), nper, object@neventsubTrt0, object@neventsubTrt1), ncol = 7)
write(" Covariate Summary Sample Type 1 Events", file = "")
write(" Subpopulation Median Minimum Maximum Size Trt Group 1 Trt Group 2", file = "")
for (i in 1:object@nsubpop) {
write(paste(format(temp[i, 1], width = 8), format(temp[i, 2], width = 15, nsmall = 2),
format(temp[i, 3], width = 12, nsmall = 4), format(temp[i, 4], width = 12, nsmall = 4),
format(temp[i, 5], width = 12), format(temp[i, 6], width = 12),
format(temp[i, 7], width = 13)), file = "")
}
} else {
write(paste("Number of subpopulations created:", object@nsubpop), file = "")
cat("\n")
write("Subpopulation summary information (including all treatments)", file = "")
nper <- apply(object@subpop, 2 ,sum)
temp <- matrix(c(1:object@nsubpop, object@medianz, round(object@minc, digits = 4),
round(object@maxc, digits = 4), nper), ncol = 5)
write(" Covariate Summary Sample", file = "")
write(" Subpopulation Median Minimum Maximum size", file = "")
for (i in 1:object@nsubpop) {
if (object@win@type == "tail-oriented" & i == length(object@win@r1) + 1 & length(object@win@r1) > 1) {
write(paste(format(temp[i, 1], width = 12), format(temp[i, 2], width = 19, nsmall = 2),
format(temp[i, 3], width = 13, nsmall = 4), format(temp[i, 4], width = 13, nsmall = 4),
format(temp[i, 5], width = 13), "(entire cohort)"), file = "")
} else if (object@win@type == "tail-oriented" & i == 1 & length(object@win@r2) > 1) {
write(paste(format(temp[i, 1], width = 12), format(temp[i, 2], width = 19, nsmall = 2),
format(temp[i, 3], width = 13, nsmall = 4), format(temp[i, 4], width = 13, nsmall = 4),
format(temp[i, 5], width = 13), "(entire cohort)"), file = "")
} else {
write(paste(format(temp[i, 1], width = 12), format(temp[i, 2], width = 19, nsmall = 2),
format(temp[i, 3], width = 13, nsmall = 4), format(temp[i, 4], width = 13, nsmall = 4),
format(temp[i, 5], width = 13)), file = "")
}
}
}
} else {
write("Subpopulations have not been generated yet.", file = "")
}
}
)
# constructor and accessor functions for stepp subpopulation
# 1. constructor
stepp.subpop <- function(swin, cov, coltype = NULL, coltrt = NULL, trts = NULL, minsubpops = NULL) {
subp <- new("stsubpop")
subp <- generate(subp, win = swin, covariate = cov, coltype = coltype, coltrt = coltrt, trts = trts,
minsubpops = minsubpops)
return(subp)
}
# 2. merge the subpopulation according to the merge vector mv
stepp.merge.subpop <- function(subpop, mv) {
subp.new <- merge(subpop, mv)
return(subp.new)
}
#3. create the edge subpopulation for the "L" or "R" side of the jth subgroup
stepp.edge.subpop <- function(subpop, j, side) {
subp.new <- edge(subpop, j, side)
return(subp.new)
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.