pred_grid_create <- function(model,data, 
                             ..., 
                             year = unique(data$year),
                             time = (unique(year) - center_time)/10,
                             center_time = mean(range(data$year))
    ) {
    time <- force(time)
    grid.args <- eval(substitute(list(...,time=time,party="a")),data)
    predigrid <- do.call(`expand.grid`,grid.args)
    
    predigrid$id <- seq.int(nrow(predigrid))
    
    predigrid %<>% within({
        CHES_id <- 1
        essround <- 3
        country <- "Germany"
    })
    for(varname in c("lrecon","religious_principles","immigration", "sociallifestyle", "galtan")) {
        if(varname %in% names(predigrid)) {
            predigrid[[paste0(varname,".sq")]] <- predigrid[[varname]]^2    
        } 
    }
    pred <- predict(model,newdata=predigrid,conditional=FALSE,se.fit=TRUE)
    
    predigrid$pred <- pred$fit
    predigrid$pred.upper <- pred$fit + 2*pred$se.fit
    predigrid$pred.lower <- pred$fit - 2*pred$se.fit

    predigrid %<>% within({
        oesch8 <- factor(oesch8, levels=levels(oesch8)[c(1,7,5,6,3,4,2,8)])
        year <- 10*time + center_time
    })
    for(varname in c("lrecon","religious_principles","immigration", "sociallifestyle", "galtan")) {
        if(varname %in% names(predigrid)) {
            predigrid[[varname]] <- 10*predigrid[[varname]] + 5
        } 
    }
    predigrid
}

prob_grid_create <- function(model,data, 
                             ..., 
                             year = unique(data$year),
                             time = (unique(year) - center_time)/10,
                             center_time = mean(range(data$year)),
                             reference = c(religious_princples=-.25),
                             transformations=NULL
    ) {
    time <- force(time)
    grid.args <- eval(substitute(list(...,time=time,party="a")),data)
    
    predigrid <- do.call(`expand.grid`,grid.args)
    transformations <- eval(substitute(transformations),predigrid)
    if(length(transformations)) {
        for(nm in names(transformations)) {
            predigrid[[nm]] <- transformations[[nm]]
        }
    }
    
    n <- nrow(predigrid)
    i <- 1:n
    predigrid$id <- i

    m <- length(reference)
    ii <- rep(i, each=2)
    predigrid <- predigrid[ii,,drop=FALSE]
    jj <- i*2 - 1
    for(vname in names(reference)) {
        predigrid[[vname]][jj] <- reference[vname]
    }
    
    predigrid %<>% within({
        CHES_id <- 1
        essround <- 3
        country <- "Germany"
    })
    for(varname in c("lrecon","religious_principles","immigration", "sociallifestyle", "galtan")) {
        if(varname %in% names(predigrid)) {
            predigrid[[paste0(varname,".sq")]] <- predigrid[[varname]]^2    
        } 
    }
    
    pred <- predict(model,newdata=predigrid,conditional=FALSE,se.fit=TRUE,type="response")
    
    predigrid$prob <- pred$fit
    predigrid$prob.upper <- pred$fit + 2*pred$se.fit
    predigrid$prob.lower <- pred$fit - 2*pred$se.fit

    predigrid %<>% within({
        oesch8 <- factor(oesch8, levels=levels(oesch8)[c(1,7,5,6,3,4,2,8)])
        year <- 10*time + center_time
    })
    for(varname in c("lrecon","religious_principles","immigration", "sociallifestyle", "galtan")) {
        if(varname %in% names(predigrid)) {
            predigrid[[varname]] <- 10*predigrid[[varname]] + 5
        } 
    }
    predigrid[i*2,,drop=FALSE]
}




ggpoints <- gtools::defmacro(data,x,y,ymin,ymax,labels=NULL,
                             percentages=FALSE,
                             ylim=NULL,
                             expr={
    if(length(ylim)) {
        ylim_ <- ylim
    }
    else {
        if(percentages){
            ylim_ <- c(0,1)
        }
        else {
            ylim_ <- c(-1.5,1.5)
        }
    }
    g <- (
        ggplot(data,aes(x=x))  
        + geom_point(aes(y=y))
        + geom_linerange(aes(ymin=ymin,ymax=ymax))
        + do.call(labs,as.list(labels))
        + coord_cartesian(ylim=ylim_) 
        + theme_bw()
        + theme(plot.title=element_text(hjust=.5),
            axis.text.x=element_text(angle=90,vjust=.5,hjust=1))) 
    if(percentages) {
        g <- g + scale_y_continuous(labels = scales::label_percent())
    }
    return(g)
})

drop_annot_y <- function(gg) {
    gg + theme(axis.title.y=element_blank(), 
               axis.text.y=element_blank(), 
               axis.ticks.y=element_blank())
}


ggcurves <- gtools::defmacro(data,x,y,ymin,ymax,grp,labels=NULL,
                             ylim=NULL,
                             percentages=FALSE,
                             vline=NULL,
                             expr={
    if(length(ylim)) {
        ylim_ <- ylim
    }
    else {
        if(percentages){
            ylim_ <- c(0,1)
        }
        else {
            ylim_ <- c(-1.5,1.5)
        }
    }
    g <- ggplot(data,aes(x=x))
    if(length(vline)) {
        g <- g + geom_vline(xintercept=vline,color="gray80")
    }
    g <- (g + geom_line(aes(y=y,linetype=grp))
            + geom_ribbon(aes(ymin=ymin,ymax=ymax,group=grp),alpha=.1)
            + do.call(labs,as.list(labels))
            + coord_cartesian(ylim=ylim_) 
            + theme_bw() + theme(legend.position="bottom")) 
    if(percentages) {
        g <- g + scale_y_continuous(labels = scales::label_percent())
    }
    return(g)
})