Skip to content

Commit

Permalink
checkpoint pre namespace changes
Browse files Browse the repository at this point in the history
  • Loading branch information
karchjd committed Mar 26, 2018
1 parent ff33559 commit ad0101e
Show file tree
Hide file tree
Showing 14 changed files with 403 additions and 91 deletions.
2 changes: 2 additions & 0 deletions .gitignore
Original file line number Diff line number Diff line change
Expand Up @@ -8,3 +8,5 @@ README.html
NAMESPACE
*.pdf
*.tex
vignettes/cache
vignettes/gppm.bbl
136 changes: 80 additions & 56 deletions R/extractors.R
Original file line number Diff line number Diff line change
@@ -1,49 +1,42 @@
#' Extract point estimates from a Gaussian process panel model
#'
#' @export
coef.GPPM <- function (object, ...){
if (is.null(object$fitRes$paraEsts)){
coef.GPPM <- function (model, ...){
if (is.null(model$fitRes$paraEsts)){
stop('Fit the model first using fit()')
}else{
object$fitRes$paraEsts
model$fitRes$paraEsts
}
}

#' @export
coef.GPPM <- function (object, ...){
if (is.null(object$fitRes$paraEsts)){
stop('Fit the model first using fit()')
}else{
object$fitRes$paraEsts
}
}

#' @export
vcov.GPPM <- function (object, ...)
vcov.GPPM <- function (model, ...)
{
object$fitRes$vcov
model$fitRes$vcov
}


#' @export
vcov.GPPM <- function (object, ...)
{
object$fitRes$vcov
SE <- function (model, ...){
UseMethod("SE")
}


#' @export
SE <- function (object, ...){
UseMethod("SE")
SE.default <- function (model, ...)
{
sqrt(diag(vcov(model, ...)))
}

#' @export
SE.default <- function (object, ...)
SE.GPPM <- function (model, ...)
{
sqrt(diag(vcov(object, ...)))
SE.default(model, ...) #here to have SE appear in extractors('GPPM') output
}

#' @export
confint.GPPM <- function (object, parm, level = 0.95, ...)
confint.GPPM <- function (model, parm, level = 0.95, ...)
{
cf <- coef(object)
cf <- coef(model)
pnames <- names(cf)
if (missing(parm))
parm <- pnames
Expand All @@ -55,24 +48,24 @@ confint.GPPM <- function (object, parm, level = 0.95, ...)
pct <- stats:::format.perc(a, 3)
ci <- array(NA, dim = c(length(parm), 2L), dimnames = list(parm,
pct))
ses <- sqrt(diag(vcov(object)))[parm]
ses <- sqrt(diag(vcov(model)))[parm]
ci[] <- cf[parm] + ses %o% fac
ci
}


#' @export
logLik.GPPM <- function(object){
val <- -2*object$fitRes$minus2LL
attr(val, "nobs") <- object$dataForStan$nPer
attr(val, "df") <- object$fitRes$nPar
logLik.GPPM <- function(model){
val <- -2*model$fitRes$minus2LL
attr(val, "nobs") <- model$dataForStan$nPer
attr(val, "df") <- model$fitRes$nPar
class(val) <- "logLik"
val
}

#' @export
fitted.GPPM <- function(object,..){
list(mean=object$fitRes$mu,cov=object$fitRes$Sigma,IDs=attr(object$dataForStan,'IDs'))
fitted.GPPM <- function(model,..){
list(mean=model$fitRes$mu,cov=model$fitRes$Sigma,IDs=attr(model$dataForStan,'IDs'))
}


Expand All @@ -82,7 +75,7 @@ nobs.GPPM <- function (x, ...) {
}

#' @export
npar <- function(theModel) {
npar <- function(model) {
UseMethod("npar")
}

Expand All @@ -92,7 +85,7 @@ npar.GPPM <- function (x, ...) {
}

#' @export
maxtime <- function(theModel) {
maxtime <- function(model) {
UseMethod("maxtime")
}

Expand All @@ -102,7 +95,7 @@ maxtime.GPPM <- function (x, ...) {
}

#' @export
nstime <- function(theModel) {
nstime <- function(model) {
UseMethod("nstime")
}

Expand All @@ -112,19 +105,19 @@ nstime.GPPM <- function (x, ...) {
}

#' @export
npred <- function(theModel) {
npred <- function(model) {
UseMethod("npred")
}

#' @export
preds <- function(theModel) {
preds <- function(model) {
UseMethod("preds")
}


#' @export
preds.GPPM <- function(theModel) {
theModel$parsedModel$preds
preds.GPPM <- function(model) {
model$parsedModel$preds
}
#' @export
npred.GPPM <- function (x, ...) {
Expand All @@ -137,32 +130,37 @@ case.names.GPPM <- function (x, ...) {
}

#' @export
variable.names.GPPM <- function (object, ...) {
object$parsedModel$params
variable.names.GPPM <- function (model, ...) {
model$parsedModel$params
}

#' @export
paramEsts <- function (model, level=.95) {
UseMethod("paramEsts")
}


#' @export
parameterEsts <- function (object, level=.95) {
paramEsts.GPPM <- function (model, level=.95) {
stopifnot(level<=1 & level>0)
stopifnot(isFitted(object))
res <- as.data.frame(matrix(nrow=object$fitRes$nPar,ncol=5))
stopifnot(isFitted(model))
res <- as.data.frame(matrix(nrow=model$fitRes$nPar,ncol=5))
a <- (1 - level)/2
a <- c(a, 1 - a)
pct <- stats:::format.perc(a, 3)
names(res) <- c('Param','Est','SE',pct)
paraNames <- variable.names(object)
confInters <- confint(object)
paraNames <- variable.names(model)
confInters <- confint(model)
res[['Param']] <- paraNames
res[['Est']] <- coef(object)[paraNames]
res[['SE']] <- SE(object)[paraNames]
res[['Est']] <- coef(model)[paraNames]
res[['SE']] <- SE(model)[paraNames]
res[[4]] <- confInters[,1][paraNames]
res[[5]] <- confInters[,2][paraNames]
res
}

#' @export
meanf <- function(theModel) {
meanf <- function(model) {
UseMethod("meanf")
}

Expand All @@ -172,7 +170,7 @@ meanf.GPPM <- function (x, ...) {
}

#' @export
covf <- function(theModel) {
covf <- function(model) {
UseMethod("covf")
}

Expand All @@ -182,23 +180,49 @@ covf.GPPM <- function (x, ...) {
}

#' @export
getIntern <- function (theModel, value, ...) {
getIntern <- function (model, value, ...) {
UseMethod("getIntern")
}

#' @export
getIntern.GPPM <- function (theModel, value, ...) {
switch(value,data=theModel$data,parsedmFormula=theModel$parsedModel$mFormula,parsedcFormula=theModel$parsedModel$kFormula,stanData=theModel$dataForStan,
stanModel=theModel$stanModel,stanOut=theModel$stanOut)
getIntern.GPPM <- function (model, value, ...) {
switch(value,data=model$data,parsedmFormula=model$parsedModel$mFormula,parsedcFormula=model$parsedModel$kFormula,stanData=model$dataForStan,
stanModel=model$stanModel,stanOut=model$stanOut)
}

#' @export
isFitted <- function(theModel) {
getData <- function (model) {
UseMethod("getData")
}

#' @export
getData.GPPM <- function (model) {
model$data
}

#' @export
getIntern.GPPM <- function (model, value, ...) {
switch(value,data=model$data,parsedmFormula=model$parsedModel$mFormula,parsedcFormula=model$parsedModel$kFormula,stanData=model$dataForStan,
stanModel=model$stanModel,stanOut=model$stanOut)
}

#' @export
isFitted <- function(model) {
UseMethod("isFitted")
}


#' @export
isFitted.GPPM <- function(theModel) {
class(theModel$fitRes)=="StanData"
isFitted.GPPM <- function(model) {
class(model$fitRes)=="StanData"
}


#' @export
extractors <- function(class) {
theMethods <- methods(class=class)
for (i in 1:length(theMethods)){

}
}

4 changes: 2 additions & 2 deletions R/fit.R
Original file line number Diff line number Diff line change
Expand Up @@ -28,7 +28,7 @@ fit <- function(theModel,...) {
#' linearChangeFake <- fit(linearChange,init=startVals,useOptimizer=FALSE)
#' stopifnot(identical(startVals,coef(linearChangeFake)))
#' @export
fit.GPPM <- function(gpModel,init='random',useOptimizer=TRUE,verbose=TRUE) {
fit.GPPM <- function(gpModel,init='random',useOptimizer=TRUE,verbose=FALSE) {
if (useOptimizer){
iter<- 10000
algorithm <- 'LBFGS' #default
Expand All @@ -44,7 +44,7 @@ fit.GPPM <- function(gpModel,init='random',useOptimizer=TRUE,verbose=TRUE) {
if(verbose){
outf <- eval
}else{
outf <- capture.output()
outf <- capture.output
}
outf(gpModel$stanOut <- rstan::optimizing(gpModel$stanModel,gpModel$dataForStan,hessian = TRUE,iter=iter,init=init,algorithm=algorithm,as_vector=FALSE))
gpModel$fitRes <- extractFitRes(gpModel$stanOut,gpModel$parsedModel,gpModel$dataForStan[c('nPer','nTime','maxTime')])
Expand Down
5 changes: 0 additions & 5 deletions R/plot.R

This file was deleted.

10 changes: 2 additions & 8 deletions R/plotLong.R
Original file line number Diff line number Diff line change
@@ -1,23 +1,17 @@
#' @import cowplot
#' @export
plotLong <- function(myData,plotIds,by,id,dv){
plot.LongData <- function(myData,plotIds,by,id,dv){
if (!is.null(attr(myData,'preds')) && missing(by)){
by <- attr(myData,'preds')
}
myData <- as_LongData(myData,id,dv)
idCol <- attr(myData,'ID')
dvCol <- attr(myData,'DV')
if(missing(plotIds)){
plotIds <- sample(unique(myData[,idCol]),3)
plotIds <- sample(unique(myData[,idCol]),5)
}
plotData <- myData[myData[,idCol] %in% plotIds,]
plotData[,idCol] <- as.factor(plotData[,idCol])
toPlot <- ggplot(plotData,aes_string(x=by,y=dvCol,colour=idCol)) + geom_line()
print(toPlot)
}


# return(ggplot(simDataLong,aes(x=Time,y=value,colour=variable)) + geom_line(size=2)
# + ylab('Value') + labs(colour = "Person")
# + theme(axis.text = element_text(size=20),axis.title = element_text(size=22),
# legend.text = element_text(size=20), legend.title= element_text(size=22)))
4 changes: 2 additions & 2 deletions R/summary.R
Original file line number Diff line number Diff line change
Expand Up @@ -47,7 +47,7 @@ summary.GPPM <- function (x, ...) {
modelSpecification <- new_ModelSpecification(meanf(x),covf(x),npar(x),variable.names(x),npred(x),preds(x))
dataStats <- new_DataStats(nobs(x),maxtime(x),nstime(x))
if (isFitted(x)){
parameterEstimates <- parameterEsts(x)
parameterEstimates <- paramEsts(x)
modelfit <- new_ModelFit(AIC(x),BIC(x),logLik(x))
}else{
parameterEstimates <- modelfit <- NA
Expand All @@ -60,7 +60,7 @@ roundForPrint <- function(x){
rf <- function(x){round(x,2)}

if(is.data.frame(x)){
sel <- vapply(tmp,is.numeric,FUN.VALUE=TRUE)
sel <- vapply(x,is.numeric,FUN.VALUE=TRUE)
x[,sel] <- rf(x[,sel])
}else{
x <- rf(x)
Expand Down
7 changes: 3 additions & 4 deletions data-raw/generateDemoLGCM.R
Original file line number Diff line number Diff line change
Expand Up @@ -4,22 +4,21 @@ timePointsPerPerson <- function(){round(rnorm(1,mean=8,sd=1))}
getTime <- function(timePoint){rnorm(1,mean=timePoint,sd=0.5)}
set.seed(249)
demoLGCM <- data.frame(matrix(nrow=1,ncol = 3))
names(demoLGCM) <- c('ID','t','x')
names(demoLGCM) <- c('ID','t','y')
counter <- 1
for (i in 1:numberPersons){
nTime <- timePointsPerPerson()
for (j in 1:nTime){
demoLGCM[counter,'ID'] <- i
demoLGCM[counter,'t'] <- getTime(j)
demoLGCM[counter,'x'] <- 1
demoLGCM[counter,'y'] <- 1
counter <- counter + 1
}
}
meanf <- 'muI+muS*t'
covf <- 'varI+covIS*(t+t#)+varS*t*t#+(t==t#)*sigma'
lgcm <- gppm(meanf,covf,demoLGCM,'ID','x')
lgcm <- gppm(meanf,covf,demoLGCM,'ID','y')
trueParas <- c(58,-1,5,1,0, 0.01)
names(trueParas) <-c('muI','muS','varI','varS','covIS','sigma')
demoLGCM <- simulate(lgcm,trueParas)
devtools::use_data(demoLGCM,trueParas,overwrite = TRUE)

Binary file modified data/demoLGCM.rda
Binary file not shown.
11 changes: 11 additions & 0 deletions man/coef.GPPM.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

3 changes: 2 additions & 1 deletion man/fit.GPPM.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

2 changes: 0 additions & 2 deletions tests/testthat/test-plotLong.R

This file was deleted.

2 changes: 1 addition & 1 deletion tests/testthat/test-simulate.R
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@ context("simulate")
test_that("working", {
simData <- simulate(lgcm,parameterValues)
yHat <- simData[,'t']*parameterValues[2]
expect_true(all.equal(yHat,simData[,'x'],tolerance=6/sqrt(nrow(simData)),scale=1))
expect_true(all.equal(yHat,simData[,3],tolerance=6/sqrt(nrow(simData)),scale=1))
})

test_that("no GPPM model as input", {
Expand Down
Loading

0 comments on commit ad0101e

Please sign in to comment.