Created
June 15, 2015 18:32
-
-
Save alexchinco/deddde7dac04f7e96c43 to your computer and use it in GitHub Desktop.
Kyle vs Grossman-Stiglitz comparative statics.
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| ## Prep workspace | |
| options(width=120, digits=6, digits.secs=6) | |
| rm(list=ls()) | |
| library(foreign) | |
| library(grid) | |
| library(plyr) | |
| library(ggplot2) | |
| library(tikzDevice) | |
| print(options('tikzLatexPackages')) | |
| options(tikzLatexPackages = | |
| c("\\usepackage{tikz}\n", | |
| "\\usepackage[active,tightpage,psfixbb]{preview}\n", | |
| "\\PreviewEnvironment{pgfpicture}\n", | |
| "\\setlength\\PreviewBorder{0pt}\n", | |
| "\\usepackage{amsmath}\n", | |
| "\\usepackage{xfrac}\n" | |
| ) | |
| ) | |
| setTikzDefaults(overwrite = FALSE) | |
| print(options('tikzLatexPackages')) | |
| library(reshape) | |
| library(foreach) | |
| ## library(doMC) | |
| ## registerDoMC(8) | |
| ## Define directories | |
| scl.str.DAT_DIR <- "~/Dropbox/research/technical_notes/" | |
| scl.str.FIG_DIR <- scl.str.DAT_DIR | |
| ## Define global baseline parameters | |
| scl.flt.SIGV <- 1 | |
| scl.flt.MUZ <- 1 | |
| scl.flt.SIGZ <- 1 | |
| scl.flt.LAM <- 1/2 | |
| scl.flt.VPHI <- 1 | |
| scl.flt.PHI <- 1 | |
| ## Define equilibrium functions | |
| fun.BET <- function(LAM, TET, PHI, SIGV) { | |
| scl.flt.OUT <- 2 * (1 - LAM) * TET * PHI * SIGV^2 | |
| return(scl.flt.OUT) | |
| } | |
| fun.KAP <- function(BET, TET, LAM, SIGV, SIGZ) { | |
| scl.flt.NUM <- TET^2 * LAM * SIGV^2 | |
| scl.flt.DEN <- TET^2 * LAM * SIGV^2 + BET^2 * SIGZ^2 | |
| scl.flt.OUT <- scl.flt.NUM/scl.flt.DEN | |
| return(scl.flt.OUT) | |
| } | |
| fun.XI <- function(PHI, VPHI, KAP, LAM, TET) { | |
| scl.flt.A <- (VPHI/PHI) * (1 - KAP * LAM) | |
| scl.flt.B <- (KAP - TET)/TET * (1 - LAM) | |
| scl.flt.OUT <- scl.flt.A - scl.flt.B | |
| return(scl.flt.OUT) | |
| } | |
| fun.TET <- function(PHI, VPHI, KAP, LAM, XI) { | |
| scl.flt.OUT <- VPHI/PHI * (1 - KAP * LAM)/XI | |
| return(scl.flt.OUT) | |
| } | |
| fun.ETA <- function(KAP, LAM, XI, TET, BET, MUZ) { | |
| scl.flt.NUM <- KAP * (1 - LAM) | |
| scl.flt.DEN <- XI * TET + KAP * (1 - LAM) | |
| scl.flt.OUT <- scl.flt.NUM/scl.flt.DEN * BET * MUZ | |
| return(scl.flt.OUT) | |
| } | |
| ## Define optimization function | |
| fun.OPTIMIZE <- function(x, vec.flt.PARAM) { | |
| LAM <- vec.flt.PARAM[1] | |
| PHI <- vec.flt.PARAM[2] | |
| VPHI <- vec.flt.PARAM[3] | |
| SIGV <- vec.flt.PARAM[4] | |
| SIGZ <- vec.flt.PARAM[5] | |
| scl.flt.BET <- fun.BET(LAM, x, PHI, SIGV) | |
| scl.flt.KAP <- fun.KAP(scl.flt.BET, x, LAM, SIGV, SIGZ) | |
| scl.flt.XI <- fun.XI(PHI, VPHI, scl.flt.KAP, LAM, x) | |
| scl.flt.TET2 <- fun.TET(PHI, VPHI, scl.flt.KAP, LAM, scl.flt.XI) | |
| scl.flt.OUT <- (x - scl.flt.TET2)^2 | |
| return(scl.flt.OUT) | |
| } | |
| fun.OPTIMIZE2 <- function(x, vec.flt.PARAM) { | |
| LAM <- vec.flt.PARAM[1] | |
| PHI <- vec.flt.PARAM[2] | |
| VPHI <- vec.flt.PARAM[3] | |
| SIGZ <- vec.flt.PARAM[4] | |
| scl.flt.SIGV <- sqrt(LAM)/(1 - LAM) * (1 - x)/(SIGZ * PHI) | |
| scl.flt.BET <- fun.BET(LAM, x, PHI, scl.flt.SIGV) | |
| scl.flt.KAP <- fun.KAP(scl.flt.BET, x, LAM, scl.flt.SIGV, SIGZ) | |
| scl.flt.XI <- fun.XI(PHI, VPHI, scl.flt.KAP, LAM, x) | |
| scl.flt.TET2 <- fun.TET(PHI, VPHI, scl.flt.KAP, LAM, scl.flt.XI) | |
| if (scl.flt.TET2 > 0) { | |
| scl.flt.OUT <- (x - scl.flt.TET2)^2 | |
| } else { | |
| scl.flt.OUT <- 10 * (x - scl.flt.TET2)^2 | |
| } | |
| return(scl.flt.OUT) | |
| } | |
| ## Sigma V | |
| if (TRUE == TRUE) { | |
| vec.flt.SIGZ <- seq(0.01, 2, by = 0.01) | |
| scl.flt.SIGZ_LEN <- length(vec.flt.SIGZ) | |
| mat.flt.PLOT <- foreach(s=1:scl.flt.SIGZ_LEN, .combine = "rbind") %do% { | |
| vec.flt.TET <- optimize(fun.OPTIMIZE2, | |
| c(0.50, 1), | |
| vec.flt.PARAM = | |
| c(scl.flt.LAM, | |
| scl.flt.PHI, | |
| scl.flt.VPHI, | |
| vec.flt.SIGZ[s] | |
| ) | |
| ) | |
| scl.flt.TET <- as.numeric(vec.flt.TET[1]) | |
| return(c(vec.flt.SIGZ[s], | |
| sqrt(scl.flt.LAM)/(1 - scl.flt.LAM) * (1 - scl.flt.TET)/(vec.flt.SIGZ[s] * scl.flt.PHI) | |
| ) | |
| ) | |
| } | |
| mat.df.PLOT <- data.frame(mat.flt.PLOT) | |
| names(mat.df.PLOT) <- c("sigz", "sigv") | |
| mat.df.PLOT <- melt(mat.df.PLOT, c("sigz")) | |
| ## Plot results | |
| theme_set(theme_bw()) | |
| scl.str.RAW_FILE <- 'plot--grossman-stiglitz-vs-kyle--tuned-sigv' | |
| scl.str.TEX_FILE <- paste(scl.str.RAW_FILE,'.tex',sep='') | |
| scl.str.PDF_FILE <- paste(scl.str.RAW_FILE,'.pdf',sep='') | |
| scl.str.PNG_FILE <- paste(scl.str.RAW_FILE,'.png',sep='') | |
| scl.str.AUX_FILE <- paste(scl.str.RAW_FILE,'.aux',sep='') | |
| scl.str.LOG_FILE <- paste(scl.str.RAW_FILE,'.log',sep='') | |
| tikz(file = scl.str.TEX_FILE, height = 2, width = 7, standAlone=TRUE) | |
| obj.gg2.PLOT <- ggplot(data = mat.df.PLOT) | |
| obj.gg2.PLOT <- obj.gg2.PLOT + scale_colour_brewer(palette = "Set1") | |
| obj.gg2.PLOT <- obj.gg2.PLOT + geom_path(aes(x = sigz, | |
| y = value * 100 | |
| ), | |
| size = 1.25 | |
| ) | |
| ## obj.gg2.PLOT <- obj.gg2.PLOT + scale_y_continuous(limits = c(0, 3), | |
| ## breaks = c(0, 1, 2, 3) | |
| ## ) | |
| obj.gg2.PLOT <- obj.gg2.PLOT + xlab("$\\sigma_z$ ${\\scriptstyle [\\text{units:}\\text{sh}]}$") | |
| obj.gg2.PLOT <- obj.gg2.PLOT + ylab("$100 \\times \\sigma_v$ ${\\scriptstyle [\\text{units:}\\sfrac{\\mathdollar 100}{\\text{sh}}]}$") | |
| obj.gg2.PLOT <- obj.gg2.PLOT + theme(plot.margin = unit(c(1,0.15,0.15,0.15), "lines"), | |
| axis.text = element_text(size = 10), | |
| axis.title = element_text(size = 10), | |
| plot.title = element_text(vjust = 1.75), | |
| panel.grid.minor = element_blank(), | |
| legend.position = "none" | |
| ) | |
| obj.gg2.PLOT <- obj.gg2.PLOT + ggtitle("Equating Demand Responses") | |
| print(obj.gg2.PLOT) | |
| dev.off() | |
| system(paste('lualatex', file.path(scl.str.TEX_FILE)), ignore.stdout = TRUE) | |
| system(paste('rm ', scl.str.TEX_FILE, sep = '')) | |
| system(paste('mv ', scl.str.PDF_FILE, ' ', scl.str.FIG_DIR, sep = '')) | |
| system(paste('rm ', scl.str.AUX_FILE, sep = '')) | |
| system(paste('rm ', scl.str.LOG_FILE, sep = '')) | |
| } | |
| ## Demand response | |
| if (TRUE == TRUE) { | |
| vec.flt.LAM <- seq(0.01, 0.99, by = 0.01) | |
| scl.flt.LAM_LEN <- length(vec.flt.LAM) | |
| mat.flt.PLOT_LAM <- foreach(l=1:scl.flt.LAM_LEN, .combine = "rbind") %do% { | |
| vec.flt.TET <- optimize(fun.OPTIMIZE, | |
| c(0, 10), | |
| vec.flt.PARAM = | |
| c(vec.flt.LAM[l], | |
| scl.flt.PHI, | |
| scl.flt.VPHI, | |
| scl.flt.SIGV, | |
| scl.flt.SIGZ | |
| ) | |
| ) | |
| scl.flt.TET <- as.numeric(vec.flt.TET[1]) | |
| return(c(vec.flt.LAM[l], | |
| 1/sqrt(vec.flt.LAM[l]) * (scl.flt.SIGZ/scl.flt.SIGV), | |
| (1/scl.flt.PHI) * (1/scl.flt.SIGV^2) * (1 - scl.flt.TET)/(1 - vec.flt.LAM[l]) | |
| ) | |
| ) | |
| } | |
| mat.df.PLOT_LAM <- data.frame(mat.flt.PLOT_LAM) | |
| names(mat.df.PLOT_LAM) <- c("lam", "Kyle", "Grossman-Stiglitz") | |
| mat.df.PLOT_LAM <- melt(mat.df.PLOT_LAM, c("lam")) | |
| vec.flt.PHI <- seq(0.01, 2, by = 0.01) | |
| scl.flt.PHI_LEN <- length(vec.flt.PHI) | |
| mat.flt.PLOT_PHI <- foreach(p=1:scl.flt.PHI_LEN, .combine = "rbind") %do% { | |
| vec.flt.TET <- optimize(fun.OPTIMIZE, | |
| c(0, 10), | |
| vec.flt.PARAM = | |
| c(scl.flt.LAM, | |
| vec.flt.PHI[p], | |
| scl.flt.VPHI, | |
| scl.flt.SIGV, | |
| scl.flt.SIGZ | |
| ) | |
| ) | |
| scl.flt.TET <- as.numeric(vec.flt.TET[1]) | |
| return(c(vec.flt.PHI[p], | |
| 1/sqrt(scl.flt.LAM) * (scl.flt.SIGZ/scl.flt.SIGV), | |
| (1/vec.flt.PHI[p]) * (1/scl.flt.SIGV^2) * (1 - scl.flt.TET)/(1 - scl.flt.LAM) | |
| ) | |
| ) | |
| } | |
| mat.df.PLOT_PHI <- data.frame(mat.flt.PLOT_PHI) | |
| names(mat.df.PLOT_PHI) <- c("phi", "Kyle", "Grossman-Stiglitz") | |
| mat.df.PLOT_PHI <- melt(mat.df.PLOT_PHI, c("phi")) | |
| vec.flt.SIGZ <- seq(0.01, 2, by = 0.01) | |
| scl.flt.SIGZ_LEN <- length(vec.flt.SIGZ) | |
| mat.flt.PLOT_SIGZ <- foreach(s=1:scl.flt.SIGZ_LEN, .combine = "rbind") %do% { | |
| vec.flt.TET <- optimize(fun.OPTIMIZE, | |
| c(0, 10), | |
| vec.flt.PARAM = | |
| c(scl.flt.LAM, | |
| scl.flt.PHI, | |
| scl.flt.VPHI, | |
| scl.flt.SIGV, | |
| vec.flt.SIGZ[s] | |
| ) | |
| ) | |
| scl.flt.TET <- as.numeric(vec.flt.TET[1]) | |
| return(c(vec.flt.SIGZ[s], | |
| 1/sqrt(scl.flt.LAM) * (vec.flt.SIGZ[s]/scl.flt.SIGV), | |
| (1/scl.flt.PHI) * (1/scl.flt.SIGV^2) * (1 - scl.flt.TET)/(1 - scl.flt.LAM) | |
| ) | |
| ) | |
| } | |
| mat.df.PLOT_SIGZ <- data.frame(mat.flt.PLOT_SIGZ) | |
| names(mat.df.PLOT_SIGZ) <- c("sigz", "Kyle", "Grossman-Stiglitz") | |
| mat.df.PLOT_SIGZ <- melt(mat.df.PLOT_SIGZ, c("sigz")) | |
| ## Plot results | |
| obj.vplayout <- function(x, y) viewport(layout.pos.row = x, layout.pos.col = y) | |
| theme_set(theme_bw()) | |
| scl.str.RAW_FILE <- 'plot--grossman-stiglitz-vs-kyle--demand-response' | |
| scl.str.TEX_FILE <- paste(scl.str.RAW_FILE,'.tex',sep='') | |
| scl.str.PDF_FILE <- paste(scl.str.RAW_FILE,'.pdf',sep='') | |
| scl.str.PNG_FILE <- paste(scl.str.RAW_FILE,'.png',sep='') | |
| scl.str.AUX_FILE <- paste(scl.str.RAW_FILE,'.aux',sep='') | |
| scl.str.LOG_FILE <- paste(scl.str.RAW_FILE,'.log',sep='') | |
| tikz(file = scl.str.TEX_FILE, height = 2, width = 7, standAlone=TRUE) | |
| ## Lambda | |
| obj.gg2.PLOT1 <- ggplot(data = mat.df.PLOT_LAM) | |
| obj.gg2.PLOT1 <- obj.gg2.PLOT1 + scale_colour_brewer(palette = "Set1") | |
| obj.gg2.PLOT1 <- obj.gg2.PLOT1 + geom_path(aes(x = lam, | |
| y = value, | |
| group = variable, | |
| colour = variable | |
| ), | |
| size = 1.25 | |
| ) | |
| obj.gg2.PLOT1 <- obj.gg2.PLOT1 + xlab("$\\lambda$ ${\\scriptstyle [\\text{units:}1]}$") | |
| obj.gg2.PLOT1 <- obj.gg2.PLOT1 + ylab("") | |
| obj.gg2.PLOT1 <- obj.gg2.PLOT1 + scale_y_continuous(limits = c(0, 3), | |
| breaks = c(0, 1, 2, 3) | |
| ) | |
| obj.gg2.PLOT1 <- obj.gg2.PLOT1 + theme(plot.margin = unit(c(0,0.15,0.15,-0.85), "lines"), | |
| axis.text = element_text(size = 10), | |
| axis.title = element_text(size = 10), | |
| plot.title = element_blank(), | |
| panel.grid.minor = element_blank(), | |
| legend.position = c(0.65, 0.80), | |
| legend.background = element_blank(), | |
| legend.title = element_blank() | |
| ) | |
| ## Phi | |
| obj.gg2.PLOT2 <- ggplot(data = mat.df.PLOT_PHI) | |
| obj.gg2.PLOT2 <- obj.gg2.PLOT2 + scale_colour_brewer(palette = "Set1") | |
| obj.gg2.PLOT2 <- obj.gg2.PLOT2 + geom_path(aes(x = phi, | |
| y = value, | |
| group = variable, | |
| colour = variable | |
| ), | |
| size = 1.25 | |
| ) | |
| obj.gg2.PLOT2 <- obj.gg2.PLOT2 + xlab("$\\phi$ ${\\scriptstyle [\\text{units:}\\sfrac{1}{\\mathdollar}]}$") | |
| obj.gg2.PLOT2 <- obj.gg2.PLOT2 + ylab("") | |
| obj.gg2.PLOT2 <- obj.gg2.PLOT2 + scale_y_continuous(limits = c(0, 3), | |
| breaks = c(0, 1, 2, 3) | |
| ) | |
| obj.gg2.PLOT2 <- obj.gg2.PLOT2 + theme(plot.margin = unit(c(0,0.15,0.15,-0.85), "lines"), | |
| axis.text = element_text(size = 10), | |
| axis.title = element_text(size = 10), | |
| plot.title = element_blank(), | |
| panel.grid.minor = element_blank(), | |
| legend.position = "none" | |
| ) | |
| ## Sigz | |
| obj.gg2.PLOT3 <- ggplot(data = mat.df.PLOT_SIGZ) | |
| obj.gg2.PLOT3 <- obj.gg2.PLOT3 + scale_colour_brewer(palette = "Set1") | |
| obj.gg2.PLOT3 <- obj.gg2.PLOT3 + geom_path(aes(x = sigz, | |
| y = value, | |
| group = variable, | |
| colour = variable | |
| ), | |
| size = 1.25 | |
| ) | |
| obj.gg2.PLOT3 <- obj.gg2.PLOT3 + scale_y_continuous(limits = c(0, 3), | |
| breaks = c(0, 1, 2, 3) | |
| ) | |
| obj.gg2.PLOT3 <- obj.gg2.PLOT3 + xlab("$\\sigma_z$ ${\\scriptstyle [\\text{units:}\\text{sh}]}$") | |
| obj.gg2.PLOT3 <- obj.gg2.PLOT3 + ylab("") | |
| obj.gg2.PLOT3 <- obj.gg2.PLOT3 + theme(plot.margin = unit(c(0,0.15,0.15,-0.85), "lines"), | |
| axis.text = element_text(size = 10), | |
| axis.title = element_text(size = 10), | |
| plot.title = element_blank(), | |
| panel.grid.minor = element_blank(), | |
| legend.position = "none" | |
| ) | |
| ## Generate output | |
| pushViewport(viewport(layout = grid.layout(5, 3))) | |
| grid.text("Demand Response: $\\mathrm{E}\\! \\left(\\sfrac{\\partial x\\!}{\\partial\\!\\hat{f}}\\right)$ ${\\scriptstyle [\\text{units:} \\sfrac{\\text{sh}^{\\scriptscriptstyle 2}\\!}{\\mathdollar}]}$", vp = viewport(layout.pos.row = 1, layout.pos.col = 1:3), vjust = 0.20, gp = gpar(fontsize=15)) | |
| print(obj.gg2.PLOT1, vp = obj.vplayout(2:5,1)) | |
| print(obj.gg2.PLOT2, vp = obj.vplayout(2:5,2)) | |
| print(obj.gg2.PLOT3, vp = obj.vplayout(2:5,3)) | |
| dev.off() | |
| system(paste('lualatex', file.path(scl.str.TEX_FILE)), ignore.stdout = TRUE) | |
| system(paste('rm ', scl.str.TEX_FILE, sep = '')) | |
| system(paste('mv ', scl.str.PDF_FILE, ' ', scl.str.FIG_DIR, sep = '')) | |
| system(paste('rm ', scl.str.AUX_FILE, sep = '')) | |
| system(paste('rm ', scl.str.LOG_FILE, sep = '')) | |
| } | |
| ## Plot price response | |
| if (TRUE == TRUE) { | |
| vec.flt.LAM <- seq(0.01, 0.99, by = 0.01) | |
| scl.flt.LAM_LEN <- length(vec.flt.LAM) | |
| mat.flt.PLOT_LAM <- foreach(l=1:scl.flt.LAM_LEN, .combine = "rbind") %do% { | |
| vec.flt.TET <- optimize(fun.OPTIMIZE, | |
| c(0, 10), | |
| vec.flt.PARAM = | |
| c(vec.flt.LAM[l], | |
| scl.flt.PHI, | |
| scl.flt.VPHI, | |
| scl.flt.SIGV, | |
| scl.flt.SIGZ | |
| ) | |
| ) | |
| scl.flt.TET <- as.numeric(vec.flt.TET[1]) | |
| return(c(vec.flt.LAM[l], 1/2, scl.flt.TET)) | |
| } | |
| mat.df.PLOT_LAM <- data.frame(mat.flt.PLOT_LAM) | |
| names(mat.df.PLOT_LAM) <- c("lam", "Kyle", "Grossman-Stiglitz") | |
| mat.df.PLOT_LAM <- melt(mat.df.PLOT_LAM, c("lam")) | |
| vec.flt.PHI <- seq(0.01, 2, by = 0.01) | |
| scl.flt.PHI_LEN <- length(vec.flt.PHI) | |
| mat.flt.PLOT_PHI <- foreach(p=1:scl.flt.PHI_LEN, .combine = "rbind") %do% { | |
| vec.flt.TET <- optimize(fun.OPTIMIZE, | |
| c(0, 10), | |
| vec.flt.PARAM = | |
| c(scl.flt.LAM, | |
| vec.flt.PHI[p], | |
| scl.flt.VPHI, | |
| scl.flt.SIGV, | |
| scl.flt.SIGZ | |
| ) | |
| ) | |
| scl.flt.TET <- as.numeric(vec.flt.TET[1]) | |
| return(c(vec.flt.PHI[p], 1/2, scl.flt.TET)) | |
| } | |
| mat.df.PLOT_PHI <- data.frame(mat.flt.PLOT_PHI) | |
| names(mat.df.PLOT_PHI) <- c("phi", "Kyle", "Grossman-Stiglitz") | |
| mat.df.PLOT_PHI <- melt(mat.df.PLOT_PHI, c("phi")) | |
| vec.flt.SIGZ <- seq(0.01, 2, by = 0.01) | |
| scl.flt.SIGZ_LEN <- length(vec.flt.SIGZ) | |
| mat.flt.PLOT_SIGZ <- foreach(s=1:scl.flt.SIGZ_LEN, .combine = "rbind") %do% { | |
| vec.flt.TET <- optimize(fun.OPTIMIZE, | |
| c(0, 10), | |
| vec.flt.PARAM = | |
| c(scl.flt.LAM, | |
| scl.flt.PHI, | |
| scl.flt.VPHI, | |
| scl.flt.SIGV, | |
| vec.flt.SIGZ[s] | |
| ) | |
| ) | |
| scl.flt.TET <- as.numeric(vec.flt.TET[1]) | |
| return(c(vec.flt.SIGZ[s], 1/2, scl.flt.TET)) | |
| } | |
| mat.df.PLOT_SIGZ <- data.frame(mat.flt.PLOT_SIGZ) | |
| names(mat.df.PLOT_SIGZ) <- c("sigz", "Kyle", "Grossman-Stiglitz") | |
| mat.df.PLOT_SIGZ <- melt(mat.df.PLOT_SIGZ, c("sigz")) | |
| ## Plot results | |
| obj.vplayout <- function(x, y) viewport(layout.pos.row = x, layout.pos.col = y) | |
| theme_set(theme_bw()) | |
| scl.str.RAW_FILE <- 'plot--grossman-stiglitz-vs-kyle--price-response' | |
| scl.str.TEX_FILE <- paste(scl.str.RAW_FILE,'.tex',sep='') | |
| scl.str.PDF_FILE <- paste(scl.str.RAW_FILE,'.pdf',sep='') | |
| scl.str.PNG_FILE <- paste(scl.str.RAW_FILE,'.png',sep='') | |
| scl.str.AUX_FILE <- paste(scl.str.RAW_FILE,'.aux',sep='') | |
| scl.str.LOG_FILE <- paste(scl.str.RAW_FILE,'.log',sep='') | |
| tikz(file = scl.str.TEX_FILE, height = 2, width = 7, standAlone=TRUE) | |
| ## Lambda | |
| obj.gg2.PLOT1 <- ggplot(data = mat.df.PLOT_LAM) | |
| obj.gg2.PLOT1 <- obj.gg2.PLOT1 + scale_colour_brewer(palette = "Set1") | |
| obj.gg2.PLOT1 <- obj.gg2.PLOT1 + geom_path(aes(x = lam, | |
| y = value, | |
| group = variable, | |
| colour = variable | |
| ), | |
| size = 1.25 | |
| ) | |
| obj.gg2.PLOT1 <- obj.gg2.PLOT1 + xlab("$\\lambda$ ${\\scriptstyle [\\text{units:}1]}$") | |
| obj.gg2.PLOT1 <- obj.gg2.PLOT1 + ylab("") | |
| obj.gg2.PLOT1 <- obj.gg2.PLOT1 + theme(plot.margin = unit(c(0,0.15,0.15,-0.85), "lines"), | |
| axis.text = element_text(size = 10), | |
| axis.title = element_text(size = 10), | |
| plot.title = element_blank(), | |
| panel.grid.minor = element_blank(), | |
| legend.position = "none" | |
| ) | |
| ## Phi | |
| obj.gg2.PLOT2 <- ggplot(data = mat.df.PLOT_PHI) | |
| obj.gg2.PLOT2 <- obj.gg2.PLOT2 + scale_colour_brewer(palette = "Set1") | |
| obj.gg2.PLOT2 <- obj.gg2.PLOT2 + geom_path(aes(x = phi, | |
| y = value, | |
| group = variable, | |
| colour = variable | |
| ), | |
| size = 1.25 | |
| ) | |
| obj.gg2.PLOT2 <- obj.gg2.PLOT2 + xlab("$\\phi$ ${\\scriptstyle [\\text{units:}\\sfrac{1}{\\mathdollar}]}$") | |
| obj.gg2.PLOT2 <- obj.gg2.PLOT2 + ylab("") | |
| obj.gg2.PLOT2 <- obj.gg2.PLOT2 + theme(plot.margin = unit(c(0,0.15,0.15,-0.85), "lines"), | |
| axis.text = element_text(size = 10), | |
| axis.title = element_text(size = 10), | |
| plot.title = element_blank(), | |
| panel.grid.minor = element_blank(), | |
| legend.position = c(0.35, 0.35), | |
| legend.background = element_blank(), | |
| legend.title = element_blank() | |
| ) | |
| ## Sigz | |
| obj.gg2.PLOT3 <- ggplot(data = mat.df.PLOT_SIGZ) | |
| obj.gg2.PLOT3 <- obj.gg2.PLOT3 + scale_colour_brewer(palette = "Set1") | |
| obj.gg2.PLOT3 <- obj.gg2.PLOT3 + geom_path(aes(x = sigz, | |
| y = value, | |
| group = variable, | |
| colour = variable | |
| ), | |
| size = 1.25 | |
| ) | |
| obj.gg2.PLOT3 <- obj.gg2.PLOT3 + xlab("$\\sigma_z$ ${\\scriptstyle [\\text{units:}\\text{sh}]}$") | |
| obj.gg2.PLOT3 <- obj.gg2.PLOT3 + ylab("") | |
| obj.gg2.PLOT3 <- obj.gg2.PLOT3 + theme(plot.margin = unit(c(0,0.15,0.15,-0.85), "lines"), | |
| axis.text = element_text(size = 10), | |
| axis.title = element_text(size = 10), | |
| plot.title = element_blank(), | |
| panel.grid.minor = element_blank(), | |
| legend.position = "none" | |
| ) | |
| ## Generate output | |
| pushViewport(viewport(layout = grid.layout(5, 3))) | |
| grid.text("Price Response: $\\mathrm{E}\\left(\\sfrac{\\partial p\\!}{\\partial\\!\\hat{f}}\\right)$ ${\\scriptstyle [\\text{units:} 1]}$", vp = viewport(layout.pos.row = 1, layout.pos.col = 1:3), vjust = 0.20, gp = gpar(fontsize=15)) | |
| print(obj.gg2.PLOT1, vp = obj.vplayout(2:5,1)) | |
| print(obj.gg2.PLOT2, vp = obj.vplayout(2:5,2)) | |
| print(obj.gg2.PLOT3, vp = obj.vplayout(2:5,3)) | |
| dev.off() | |
| system(paste('lualatex', file.path(scl.str.TEX_FILE)), ignore.stdout = TRUE) | |
| system(paste('rm ', scl.str.TEX_FILE, sep = '')) | |
| system(paste('mv ', scl.str.PDF_FILE, ' ', scl.str.FIG_DIR, sep = '')) | |
| system(paste('rm ', scl.str.AUX_FILE, sep = '')) | |
| system(paste('rm ', scl.str.LOG_FILE, sep = '')) | |
| } | |
| ## Plot average price | |
| if (TRUE == TRUE) { | |
| vec.flt.LAM <- seq(0.01, 0.99, by = 0.01) | |
| scl.flt.LAM_LEN <- length(vec.flt.LAM) | |
| mat.flt.PLOT_LAM <- foreach(l=1:scl.flt.LAM_LEN, .combine = "rbind") %do% { | |
| vec.flt.TET <- optimize(fun.OPTIMIZE, | |
| c(0, 10), | |
| vec.flt.PARAM = | |
| c(vec.flt.LAM[l], | |
| scl.flt.PHI, | |
| scl.flt.VPHI, | |
| scl.flt.SIGV, | |
| scl.flt.SIGZ | |
| ) | |
| ) | |
| scl.flt.TET <- as.numeric(vec.flt.TET[1]) | |
| scl.flt.BET <- fun.BET(vec.flt.LAM[l], scl.flt.TET, scl.flt.PHI, scl.flt.SIGV) | |
| scl.flt.KAP <- fun.KAP(scl.flt.BET, scl.flt.TET, vec.flt.LAM[l], scl.flt.SIGV, scl.flt.SIGZ) | |
| scl.flt.XI <- fun.XI(scl.flt.PHI, scl.flt.VPHI, scl.flt.KAP, vec.flt.LAM[l], scl.flt.TET) | |
| scl.flt.ETA <- fun.ETA(scl.flt.KAP, vec.flt.LAM[l], scl.flt.XI, scl.flt.TET, scl.flt.BET, scl.flt.MUZ) | |
| return(c(vec.flt.LAM[l], 0, scl.flt.ETA - scl.flt.BET * scl.flt.MUZ)) | |
| } | |
| mat.df.PLOT_LAM <- data.frame(mat.flt.PLOT_LAM) | |
| names(mat.df.PLOT_LAM) <- c("lam", "Kyle", "Grossman-Stiglitz") | |
| mat.df.PLOT_LAM <- melt(mat.df.PLOT_LAM, c("lam")) | |
| vec.flt.PHI <- seq(0.01, 2, by = 0.01) | |
| scl.flt.PHI_LEN <- length(vec.flt.PHI) | |
| mat.flt.PLOT_PHI <- foreach(p=1:scl.flt.PHI_LEN, .combine = "rbind") %do% { | |
| vec.flt.TET <- optimize(fun.OPTIMIZE, | |
| c(0, 10), | |
| vec.flt.PARAM = | |
| c(scl.flt.LAM, | |
| vec.flt.PHI[p], | |
| scl.flt.VPHI, | |
| scl.flt.SIGV, | |
| scl.flt.SIGZ | |
| ) | |
| ) | |
| scl.flt.TET <- as.numeric(vec.flt.TET[1]) | |
| scl.flt.BET <- fun.BET(scl.flt.LAM, scl.flt.TET, vec.flt.PHI[p], scl.flt.SIGV) | |
| scl.flt.KAP <- fun.KAP(scl.flt.BET, scl.flt.TET, scl.flt.LAM, scl.flt.SIGV, scl.flt.SIGZ) | |
| scl.flt.XI <- fun.XI(vec.flt.PHI[p], scl.flt.VPHI, scl.flt.KAP, scl.flt.LAM, scl.flt.TET) | |
| scl.flt.ETA <- fun.ETA(scl.flt.KAP, scl.flt.LAM, scl.flt.XI, scl.flt.TET, scl.flt.BET, scl.flt.MUZ) | |
| return(c(vec.flt.PHI[p], 0, scl.flt.ETA - scl.flt.BET * scl.flt.MUZ)) | |
| } | |
| mat.df.PLOT_PHI <- data.frame(mat.flt.PLOT_PHI) | |
| names(mat.df.PLOT_PHI) <- c("phi", "Kyle", "Grossman-Stiglitz") | |
| mat.df.PLOT_PHI <- melt(mat.df.PLOT_PHI, c("phi")) | |
| vec.flt.SIGZ <- seq(0.01, 2, by = 0.01) | |
| scl.flt.SIGZ_LEN <- length(vec.flt.SIGZ) | |
| mat.flt.PLOT_SIGZ <- foreach(s=1:scl.flt.SIGZ_LEN, .combine = "rbind") %do% { | |
| vec.flt.TET <- optimize(fun.OPTIMIZE, | |
| c(0, 10), | |
| vec.flt.PARAM = | |
| c(scl.flt.LAM, | |
| scl.flt.PHI, | |
| scl.flt.VPHI, | |
| scl.flt.SIGV, | |
| vec.flt.SIGZ[s] | |
| ) | |
| ) | |
| scl.flt.TET <- as.numeric(vec.flt.TET[1]) | |
| scl.flt.BET <- fun.BET(scl.flt.LAM, scl.flt.TET, scl.flt.PHI, scl.flt.SIGV) | |
| scl.flt.KAP <- fun.KAP(scl.flt.BET, scl.flt.TET, scl.flt.LAM, scl.flt.SIGV, vec.flt.SIGZ[s]) | |
| scl.flt.XI <- fun.XI(scl.flt.PHI, scl.flt.VPHI, scl.flt.KAP, scl.flt.LAM, scl.flt.TET) | |
| scl.flt.ETA <- fun.ETA(scl.flt.KAP, scl.flt.LAM, scl.flt.XI, scl.flt.TET, scl.flt.BET, scl.flt.MUZ) | |
| return(c(vec.flt.SIGZ[s], 0, scl.flt.ETA - scl.flt.BET * scl.flt.MUZ)) | |
| } | |
| mat.df.PLOT_SIGZ <- data.frame(mat.flt.PLOT_SIGZ) | |
| names(mat.df.PLOT_SIGZ) <- c("sigz", "Kyle", "Grossman-Stiglitz") | |
| mat.df.PLOT_SIGZ <- melt(mat.df.PLOT_SIGZ, c("sigz")) | |
| ## Plot results | |
| obj.vplayout <- function(x, y) viewport(layout.pos.row = x, layout.pos.col = y) | |
| theme_set(theme_bw()) | |
| scl.str.RAW_FILE <- 'plot--grossman-stiglitz-vs-kyle--average-price' | |
| scl.str.TEX_FILE <- paste(scl.str.RAW_FILE,'.tex',sep='') | |
| scl.str.PDF_FILE <- paste(scl.str.RAW_FILE,'.pdf',sep='') | |
| scl.str.PNG_FILE <- paste(scl.str.RAW_FILE,'.png',sep='') | |
| scl.str.AUX_FILE <- paste(scl.str.RAW_FILE,'.aux',sep='') | |
| scl.str.LOG_FILE <- paste(scl.str.RAW_FILE,'.log',sep='') | |
| tikz(file = scl.str.TEX_FILE, height = 2, width = 7, standAlone=TRUE) | |
| ## Lambda | |
| obj.gg2.PLOT1 <- ggplot(data = mat.df.PLOT_LAM) | |
| obj.gg2.PLOT1 <- obj.gg2.PLOT1 + scale_colour_brewer(palette = "Set1") | |
| obj.gg2.PLOT1 <- obj.gg2.PLOT1 + geom_path(aes(x = lam, | |
| y = value, | |
| group = variable, | |
| colour = variable | |
| ), | |
| size = 1.25 | |
| ) | |
| obj.gg2.PLOT1 <- obj.gg2.PLOT1 + xlab("$\\lambda$ ${\\scriptstyle [\\text{units:}1]}$") | |
| obj.gg2.PLOT1 <- obj.gg2.PLOT1 + ylab("") | |
| obj.gg2.PLOT1 <- obj.gg2.PLOT1 + scale_y_continuous(limits = c(-1,0), | |
| breaks = c(-1, -0.8, -0.60, -0.40, -0.20, 0) | |
| ) | |
| obj.gg2.PLOT1 <- obj.gg2.PLOT1 + theme(plot.margin = unit(c(0,0.15,0.15,-0.85), "lines"), | |
| axis.text = element_text(size = 10), | |
| axis.title = element_text(size = 10), | |
| plot.title = element_blank(), | |
| panel.grid.minor = element_blank(), | |
| legend.position = "none" | |
| ) | |
| ## Phi | |
| obj.gg2.PLOT2 <- ggplot(data = mat.df.PLOT_PHI) | |
| obj.gg2.PLOT2 <- obj.gg2.PLOT2 + scale_colour_brewer(palette = "Set1") | |
| obj.gg2.PLOT2 <- obj.gg2.PLOT2 + geom_path(aes(x = phi, | |
| y = value, | |
| group = variable, | |
| colour = variable | |
| ), | |
| size = 1.25 | |
| ) | |
| obj.gg2.PLOT2 <- obj.gg2.PLOT2 + xlab("$\\phi$ ${\\scriptstyle [\\text{units:}\\sfrac{1}{\\mathdollar}]}$") | |
| obj.gg2.PLOT2 <- obj.gg2.PLOT2 + ylab("") | |
| obj.gg2.PLOT2 <- obj.gg2.PLOT2 + scale_y_continuous(limits = c(-1,0), | |
| breaks = c(-1, -0.8, -0.60, -0.40, -0.20, 0) | |
| ) | |
| obj.gg2.PLOT2 <- obj.gg2.PLOT2 + theme(plot.margin = unit(c(0,0.15,0.15,-0.85), "lines"), | |
| axis.text = element_text(size = 10), | |
| axis.title = element_text(size = 10), | |
| plot.title = element_blank(), | |
| panel.grid.minor = element_blank(), | |
| legend.position = c(0.35, 0.225), | |
| legend.background = element_blank(), | |
| legend.title = element_blank() | |
| ) | |
| ## Sigz | |
| obj.gg2.PLOT3 <- ggplot(data = mat.df.PLOT_SIGZ) | |
| obj.gg2.PLOT3 <- obj.gg2.PLOT3 + scale_colour_brewer(palette = "Set1") | |
| obj.gg2.PLOT3 <- obj.gg2.PLOT3 + geom_path(aes(x = sigz, | |
| y = value, | |
| group = variable, | |
| colour = variable | |
| ), | |
| size = 1.25 | |
| ) | |
| obj.gg2.PLOT3 <- obj.gg2.PLOT3 + xlab("$\\sigma_z$ ${\\scriptstyle [\\text{units:}\\text{sh}]}$") | |
| obj.gg2.PLOT3 <- obj.gg2.PLOT3 + ylab("") | |
| obj.gg2.PLOT3 <- obj.gg2.PLOT3 + scale_y_continuous(limits = c(-1,0), | |
| breaks = c(-1, -0.8, -0.60, -0.40, -0.20, 0) | |
| ) | |
| obj.gg2.PLOT3 <- obj.gg2.PLOT3 + theme(plot.margin = unit(c(0,0.15,0.15,-0.85), "lines"), | |
| axis.text = element_text(size = 10), | |
| axis.title = element_text(size = 10), | |
| plot.title = element_blank(), | |
| panel.grid.minor = element_blank(), | |
| legend.position = "none" | |
| ) | |
| ## Generate output | |
| pushViewport(viewport(layout = grid.layout(5, 3))) | |
| grid.text("Average Price: $\\mathrm{E}(p)$ ${\\scriptstyle [\\text{units:} \\sfrac{\\mathdollar}{\\text{sh}}]}$", vp = viewport(layout.pos.row = 1, layout.pos.col = 1:3), vjust = 0.20, gp = gpar(fontsize=15)) | |
| print(obj.gg2.PLOT1, vp = obj.vplayout(2:5,1)) | |
| print(obj.gg2.PLOT2, vp = obj.vplayout(2:5,2)) | |
| print(obj.gg2.PLOT3, vp = obj.vplayout(2:5,3)) | |
| dev.off() | |
| system(paste('lualatex', file.path(scl.str.TEX_FILE)), ignore.stdout = TRUE) | |
| system(paste('rm ', scl.str.TEX_FILE, sep = '')) | |
| system(paste('mv ', scl.str.PDF_FILE, ' ', scl.str.FIG_DIR, sep = '')) | |
| system(paste('rm ', scl.str.AUX_FILE, sep = '')) | |
| system(paste('rm ', scl.str.LOG_FILE, sep = '')) | |
| } | |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment