Gaussian Process Regression demonstrated

[1] "Created: Wed Sep  3 16:20:17 2014"

Here we demonstrate how Figure 2 is constructed:

Load in the data and helper functions:

load("lcdb.rda")
source("funcs.R")

Set up the start date, range of dates and a grid for computing the estimates:

daterange <- 2764
firstdate <- 53464
ngrid <- 100
grid <- seq(0,daterange,length=ngrid)

Plot the data for the first example

lcobj <- "CSS071216:110407-045134"
lc <- lcdb[[which(lcobj == lcdbi$id)]]
lc$mjd <- lc$mjd-firstdate
p1 <- ggplot(data=lc,aes(x=mjd,y=mag)) + xlab("Day")+ylab("Magnitude")+ggtitle("AGN")
p1 <- p1+xlim(0,daterange)+ylim(21.1,min(lc$mag)) + geom_point(size=1)

Set some of the prior parameters:

sig <- mean(lc$magerr)
sigv <- 0.27
wvec <- 2e-4
#kerntype <- "matern"
kerntype <- "exponential"
#wvec <- wvec * 100

Use a prior mean equal to the median, calculate and plot the posterior:

priormean <- rep(median(lc$mag), ngrid)
rmod <- gpreg(lc$mjd, lc$mag, grid, invwid=wvec, noisevar=sig^2, sigvar=sigv^2, priormean=priormean,kerntype=kerntype)
fdf <- data.frame(mjd=rmod$grid,mag=rmod$predict)
p1 <- p1 + geom_line(aes(x=mjd,y=mag),data=fdf)

Prior mean now set at the detection limit:

priormean <- rep(21, ngrid)
rmod <- gpreg(lc$mjd, lc$mag, grid, invwid=wvec, noisevar=sig^2, sigvar=sigv^2, priormean=priormean,kerntype=kerntype)
fdf <- data.frame(mjd=rmod$grid,mag=rmod$predict)
p1 <- p1 + geom_line(aes(mjd,mag),fdf,linetype=2)

Giving this final resulting plot:

p1

Repeat for the second example:

lcobj <- "CSS110405:141104+011153"
lc <- lcdb[[which(lcobj == lcdbi$id)]]
lc$mjd <- lc$mjd-firstdate
p2 <- ggplot(data=lc,aes(x=mjd,y=mag)) + xlab("Day")+ylab("Magnitude")+ggtitle("SN")
p2 <- p2+xlim(0,daterange)+ylim(21.1,min(lc$mag)) + geom_point(size=1)

sig <- mean(lc$magerr)
priormean <- rep(median(lc$mag), ngrid)
rmod <- gpreg(lc$mjd, lc$mag, grid, invwid=wvec, noisevar=sig^2, sigvar=sigv^2, priormean=priormean,kerntype=kerntype)
fdf <- data.frame(mjd=rmod$grid,mag=rmod$predict)
p2 <- p2 + geom_line(aes(x=mjd,y=mag),data=fdf)

priormean <- rep(21, ngrid)
rmod <- gpreg(lc$mjd, lc$mag, grid, invwid=wvec, noisevar=sig^2, sigvar=sigv^2, priormean=priormean,kerntype=kerntype)
fdf <- data.frame(mjd=rmod$grid,mag=rmod$predict)
p2 <- p2 + geom_line(aes(mjd,mag),fdf,linetype=2)
p2

Third example:

lcobj <- "CSS111103:230309+400608"
lc <- lcdb[[which(lcobj == lcdbi$id)]]
lc$mjd <- lc$mjd-firstdate
p3 <- ggplot(data=lc,aes(x=mjd,y=mag)) + xlab("Day")+ylab("Magnitude")+ggtitle("Flare")
p3 <- p3+xlim(0,daterange)+ylim(21.1,min(lc$mag)) + geom_point(size=1)

sig <- mean(lc$magerr)
priormean <- rep(median(lc$mag), ngrid)
rmod <- gpreg(lc$mjd, lc$mag, grid, invwid=wvec, noisevar=sig^2, sigvar=sigv^2, priormean=priormean,kerntype=kerntype)
fdf <- data.frame(mjd=rmod$grid,mag=rmod$predict)
p3 <- p3 + geom_line(aes(x=mjd,y=mag),data=fdf)

priormean <- rep(21, ngrid)
rmod <- gpreg(lc$mjd, lc$mag, grid, invwid=wvec, noisevar=sig^2, sigvar=sigv^2, priormean=priormean,kerntype=kerntype)
fdf <- data.frame(mjd=rmod$grid,mag=rmod$predict)
p3 <- p3 + geom_line(aes(mjd,mag),fdf,linetype=2)
p3

Fourth example:

lcobj <- "3019048007678"
lc <- lcdb[[which(lcobj == lcdbi$id)]]
lc$mjd <- lc$mjd-firstdate
p4 <- ggplot(data=lc,aes(x=mjd,y=mag)) + xlab("Day")+ylab("Magnitude")+ggtitle("Non Transient")
p4 <- p4+xlim(0,daterange)+ylim(21.1,min(lc$mag)) + geom_point(size=1)

sig <- mean(lc$magerr)
priormean <- rep(median(lc$mag), ngrid)
rmod <- gpreg(lc$mjd, lc$mag, grid, invwid=wvec, noisevar=sig^2, sigvar=sigv^2, priormean=priormean,kerntype=kerntype)
fdf <- data.frame(mjd=rmod$grid,mag=rmod$predict)
p4 <- p4 + geom_line(aes(x=mjd,y=mag),data=fdf)

priormean <- rep(21, ngrid)
rmod <- gpreg(lc$mjd, lc$mag, grid, invwid=wvec, noisevar=sig^2, sigvar=sigv^2, priormean=priormean,kerntype=kerntype)
fdf <- data.frame(mjd=rmod$grid,mag=rmod$predict)
p4 <- p4 + geom_line(aes(mjd,mag),fdf,linetype=2)
p4