#| '!! shinylive warning !!': |
#| shinylive does not work in self-contained HTML documents.
#| Please set `embed-resources: false` in your metadata.
#| standalone: true
#| viewerHeight: 700
#| components: [viewer]
## file: app.R
# ---------------------------------------------------------------------------
# Slicing a bivariate normal.
#
# (X, Y) is bivariate normal with means mu_X, mu_Y, SDs sigma_X, sigma_Y and
# correlation rho, all on sliders. The 3D panel draws the joint density, with
# the two marginals kept faint on the back walls so the surface stays the
# focus. A vertical slice at X = x cuts the surface, and the
# cross-section is the joint density along that line. Rescaled to area 1, that
# cross-section is the conditional density of Y, which the lower panel sets
# against Y's marginal:
#
# Y | X = x ~ N( mu_Y + rho * sigma_Y / sigma_X * (x - mu_X),
# sigma_Y^2 * (1 - rho^2) )
#
# so the conditional SD is sqrt(1 - rho^2) times the marginal SD, whatever x is.
#
# Deliberately base R + shiny only. Every extra package is another wasm
# download when this runs in the browser through shinylive. The 3D panel is
# drawn by plotly.js loaded straight from a CDN rather than the plotly R
# package (which would pull in ggplot2, data.table, httr, ...): R computes
# every number and sends the arrays over; the JS only draws them.
# ---------------------------------------------------------------------------
library(shiny)
# The 3D box is kept tight around the density so it isn't mostly flat floor;
# the mean and SD sliders are limited to what fits in it. The 2D panel has
# room to show the tails, so it runs wider.
LIM <- 3 # 3D x and y run over [-LIM, LIM]
G <- seq(-LIM, LIM, length.out = 61) # 0.1 steps, fine enough for the rho = 0.95 ridge
LIM_2D <- 5
COLS <- c(X = "#c5050c", Y = "#0479A8")
COND_COL <- "#C77400"
INK <- "#333333"
PLOTLY_JS <- "https://cdn.jsdelivr.net/npm/plotly.js-gl3d-dist-min@3.7.0/plotly-gl3d.min.js"
dbvn <- function(x, y, p) {
zx <- (x - p$mx) / p$sx
zy <- (y - p$my) / p$sy
q <- 1 - p$rho^2
exp(-(zx^2 - 2 * p$rho * zx * zy + zy^2) / (2 * q)) / (2 * pi * p$sx * p$sy * sqrt(q))
}
# Y | X = at, and the conditional-mean line E[Y | X = x] = my + slope * (x - mx)
slice_of <- function(p, at) {
slope <- p$rho * p$sy / p$sx
list(at = at, slope = slope, m_marg = p$my, s_marg = p$sy,
m_cond = p$my + slope * (at - p$mx), s_cond = p$sy * sqrt(1 - p$rho^2))
}
css <- "
@import url('https://fonts.googleapis.com/css2?family=Red+Hat+Display:wght@400;700&family=Red+Hat+Text:wght@400;500;700&display=swap');
html, body { height:100%; }
body { font-family:'Red Hat Text',system-ui,sans-serif; background:#F7F7F7; color:#333;
margin:0; padding:10px; font-size:14px; box-sizing:border-box; overflow:hidden; }
h4 { font-family:'Red Hat Display',system-ui,sans-serif; margin:0; font-size:17px; }
h5 { font-family:'Red Hat Display',system-ui,sans-serif; margin:2px 0 4px 0; font-size:14px;
font-weight:700; color:#666; }
/* fluidPage() wraps the UI in an auto-height .container-fluid; it needs a
definite height for the panels below to size to the frame. */
body > .container-fluid { height:100%; padding:0; }
.wrap { display:flex; gap:12px; height:100%; }
.card { display:flex; flex-direction:column; min-width:0; min-height:0;
background:#fff; border:1px solid #e2e2e2; border-radius:8px; padding:10px 12px;
box-sizing:border-box; }
.side { width:232px; flex:none; overflow-y:auto; }
.side hr { margin:8px 0 6px 0; border-color:#eee; }
.main { flex:1 1 0; min-width:0; display:flex; flex-direction:column; gap:12px; }
.c3d { flex:1.75 1 0; }
.c2d { flex:1 1 0; }
.hdr { display:flex; justify-content:space-between; align-items:baseline; height:24px; flex:none; }
.hdr .opts { display:flex; align-items:baseline; gap:14px; font-size:12.5px; color:#666; }
.hdr .form-group { width:auto !important; margin:0; }
.hdr .checkbox { margin:0; }
.hdr a { color:#0479A8; cursor:pointer; }
.model { font-size:12.5px; color:#666; }
.hdr h4 .m, .fit .m { color:#0479A8; }
.hdr h4 .c, .fit .c { color:#C77400; }
#bvn3d { flex:1 1 auto; min-height:0; }
#bvn3d .nojs { color:#999; margin:40px; text-align:center; }
.plot2d { flex:1 1 auto; min-height:0; }
.form-group { margin-bottom:0; }
.control-label { font-weight:500; white-space:nowrap; font-size:13px; margin-bottom:0; }
.irs--shiny .irs-bar { border-color:var(--gc); background:var(--gc); }
.irs--shiny .irs-single { background:var(--gc); }
.irs--shiny .irs-handle { border-color:var(--gc); }
.irs--shiny .irs-min, .irs--shiny .irs-max { display:none !important; }
.slider-animate-container { margin-top:-4px; font-size:13px; }
.fit { display:flex; gap:22px; font-size:12px; color:#777; margin:4px 0 0 0; height:28px;
flex:none; align-items:center; }
.fit b { font-family:'Red Hat Display',sans-serif; font-size:16px; color:#333; }
.meter { display:inline-block; width:70px; height:8px; background:#eee; border-radius:3px;
overflow:hidden; vertical-align:middle; margin-left:4px; }
.meter div { height:8px; }
"
# Draws the 3D panel from the arrays the server sends as a `bvn` message.
js <- r"---(
(function () {
var LIM = 3, WALL = LIM - 0.02; // curves sit just off the walls so they don't z-fight
var COL = { X: '#c5050c', Y: '#0479A8', slice: '#C77400', cef: '#333333' };
var FLAT = { ambient: 1, diffuse: 0, specular: 0, fresnel: 0 };
// The marginals are backdrop, thin and faint, so the joint surface and the
// slice carry the picture
var MARG = { width: 3.5, line: 0.45, fill: 0.06, label: 0.6 };
// The home view faces the slice from the +x side, so the cross-section and
// Y's marginal on the back wall are both seen head-on
var HOME = { az: -30, el: 26 };
var EYE_R = 1.25, ZOOM = 1.45, Z_ASPECT = 0.7;
var last = null, reset = false, near = null, hooked = false;
function arr(v) { return [].concat(v); }
function rep(v, n) { var a = []; for (var i = 0; i < n; i++) a.push(v); return a; }
// Shaded area between the floor and a curve, in the curve's vertical plane
function area(xs, ys, zs, color, opacity, visible) {
var n = xs.length, X = [], Y = [], Z = [], I = [], J = [], K = [];
for (var k = 0; k < n; k++) { X.push(xs[k]); Y.push(ys[k]); Z.push(0); }
for (var k = 0; k < n; k++) { X.push(xs[k]); Y.push(ys[k]); Z.push(zs[k]); }
for (var k = 0; k < n - 1; k++) { I.push(k, k + 1); J.push(k + 1, n + k + 1); K.push(n + k, n + k); }
return { type: 'mesh3d', x: X, y: Y, z: Z, i: I, j: J, k: K, color: color, opacity: opacity,
flatshading: true, lighting: FLAT, hoverinfo: 'skip', visible: visible };
}
function curve(xs, ys, zs, color, width, visible, dash, opacity) {
return { type: 'scatter3d', mode: 'lines', x: xs, y: ys, z: zs, visible: visible,
line: { color: color, width: width, dash: dash || 'solid' }, hoverinfo: 'skip',
opacity: opacity === undefined ? 1 : opacity };
}
function axis(title, range) {
return { title: { text: title }, range: range, gridcolor: '#dddddd', zeroline: false,
showspikes: false, showbackground: true, backgroundcolor: '#F4F4F4',
tickfont: { size: 11 }, tickangle: 0, nticks: 6 };
}
// A marginal's label, about one SD out from its peak. Y's goes on the low-y
// side, which from the home view is away from the surface.
function label(g, f, onX, color, text) {
var i = f.indexOf(Math.max.apply(null, f)), peak = f[i], sd = 0.3989 / peak;
var at = onX ? Math.max(g[i] - 1.1 * sd, -LIM + 0.6) : Math.min(g[i] + 1.1 * sd, LIM - 0.6);
return { x: onX ? -WALL : at, y: onX ? at : WALL, z: 0.8 * peak, text: text,
showarrow: false, xanchor: onX ? 'right' : 'left', opacity: MARG.label,
font: { size: 14, color: color } };
}
function homeCamera() {
var a = HOME.az * Math.PI / 180, e = HOME.el * Math.PI / 180;
return { eye: { x: EYE_R * Math.cos(e) * Math.cos(a), y: EYE_R * Math.cos(e) * Math.sin(a),
z: EYE_R * Math.sin(e) },
up: { x: 0, y: 0, z: 1 }, center: { x: 0, y: 0, z: -0.05 },
projection: { type: 'orthographic' } };
}
// Which side of the slice the camera is on (true = the high-x side)
function nearSide(camera) { return camera.eye.x > 0; }
// Split the surface at the slice into the part behind it, drawn solid, and
// the part between the slice and the camera, ghosted so the cut face shows.
// The grid and the slice slider share a 0.1 step, so the cut falls on a
// grid line, which both halves keep. Rows of z are y, columns x.
function cut(z, g, at, nearHigh) {
var solid = [], ghost = [];
for (var i = 0; i < g.length; i++) {
var zs = [], zg = [];
for (var j = 0; j < g.length; j++) {
var d = g[j] - at;
var edge = Math.abs(d) < 1e-6, inFront = !edge && (nearHigh ? d > 0 : d < 0);
zs.push(inFront ? null : z[i][j]);
zg.push(inFront || edge ? z[i][j] : null);
}
solid.push(zs); ghost.push(zg);
}
return { solid: solid, ghost: ghost };
}
function render(m) {
var el = document.getElementById('bvn3d');
if (!window.Plotly) {
el.innerHTML = '<p class="nojs">The 3D view needs plotly.js, which loads from the internet.</p>';
return;
}
last = m;
var g = arr(m.grid), n = g.length, s = m.slice, c = m.cef, zmax = m.zmax;
var margX = arr(m.margX), margY = arr(m.margY), sz = arr(s.z);
var slx = rep(s.at, n);
// Keep wherever the viewer has rotated / zoomed to. Read the live scene:
// plotly writes a drag back to the layout only after a delay, so a slider
// moved straight after a drag would otherwise snap the camera back.
var sc = el._fullLayout && el._fullLayout.scene && el._fullLayout.scene._scene;
var camera = homeCamera(), aspect = { x: ZOOM, y: ZOOM, z: Z_ASPECT * ZOOM };
if (sc && !reset) {
camera = sc.getCamera();
camera.projection = { type: 'orthographic' };
aspect = sc.glplot.getAspectratio();
}
reset = false;
near = nearSide(camera);
var z = s.show ? cut(m.z, g, s.at, near) : { solid: m.z, ghost: [[null]] };
// The colour scale spans the surface's own height, so its peak is always
// the darkest purple on screen. Purple: the red of X mixed with the blue of Y
var SURF = { type: 'surface', x: g, y: g, cmin: 0, cmax: m.peak, showscale: false,
hoverinfo: 'skip',
colorscale: [[0, '#e9dff1'], [0.35, '#9b6cc2'], [1, '#4a1766']],
lighting: { ambient: 0.7, diffuse: 0.5, specular: 0.1, roughness: 0.8, fresnel: 0.05 } };
var traces = [
Object.assign({}, SURF, { z: z.solid,
contours: { x: { highlight: false }, y: { highlight: false },
z: { show: true, start: zmax / 9, end: zmax, size: zmax / 9,
color: '#ffffff', width: 1, highlight: false } } }),
Object.assign({}, SURF, { z: z.ghost, opacity: 0.2, visible: s.show,
contours: { x: { highlight: false }, y: { highlight: false }, z: { highlight: false } } }),
area(g, rep(WALL, n), margX, COL.X, MARG.fill, true),
curve(g, rep(WALL, n), margX, COL.X, MARG.width, true, 'solid', MARG.line),
area(rep(-WALL, n), g, margY, COL.Y, MARG.fill, true),
curve(rep(-WALL, n), g, margY, COL.Y, MARG.width, true, 'solid', MARG.line),
// the cutting plane, as two triangles
{ type: 'mesh3d', x: rep(s.at, 4), y: [-LIM, LIM, LIM, -LIM], z: [0, 0, zmax, zmax],
i: [0, 0], j: [1, 2], k: [2, 3],
color: '#777777', opacity: 0.1, flatshading: true, lighting: FLAT, hoverinfo: 'skip',
visible: s.show },
area(slx, g, sz, COL.slice, 0.5, s.show),
curve(slx, g, sz, COL.slice, 9, s.show),
curve(arr(c.x), arr(c.y), rep(0, arr(c.x).length), COL.cef, 5, c.show, 'dash'),
{ type: 'scatter3d', mode: 'markers', x: [s.at], y: [s.mean],
z: [0], marker: { size: 5, color: COL.slice, symbol: 'diamond' }, hoverinfo: 'skip',
visible: s.show }
];
var layout = {
paper_bgcolor: '#ffffff', margin: { l: 0, r: 0, t: 0, b: 0 }, showlegend: false,
font: { family: "'Red Hat Text', sans-serif", color: '#333333' },
scene: {
xaxis: axis('x', [-LIM, LIM]), yaxis: axis('y', [-LIM, LIM]), zaxis: axis('', [0, zmax]),
aspectmode: 'manual', aspectratio: aspect, camera: camera,
dragmode: 'turntable', hovermode: false,
annotations: [label(g, margX, false, COL.X, 'f<sub>X</sub>(x)'),
label(g, margY, true, COL.Y, 'f<sub>Y</sub>(y)')]
}
};
Plotly.react(el, traces, layout, { displayModeBar: false, responsive: true, doubleClick: false });
// Rotating round to the other side of the slice moves the ghosted half
if (!hooked) {
hooked = true;
el.on('plotly_relayout', function () {
var sc = el._fullLayout.scene._scene;
if (last && nearSide(sc.getCamera()) !== near) render(last);
});
}
}
Shiny.addCustomMessageHandler('bvn', render);
document.addEventListener('click', function (e) {
if (e.target && e.target.id === 'reset3d') {
e.preventDefault();
reset = true;
if (last) render(last);
}
});
})();
)---"
slider <- function(id, label, min, max, value, step, col, ...) {
div(style = paste0("--gc:", col, ";"),
sliderInput(id, HTML(label), min = min, max = max, value = value, step = step,
width = "100%", ticks = FALSE, ...))
}
ui <- fluidPage(
tags$head(tags$style(HTML(css)), tags$script(src = PLOTLY_JS)),
div(class = "wrap",
div(class = "card side",
h5("Joint distribution"),
slider("rho", "ρ – correlation", -0.95, 0.95, 0, 0.05, INK),
slider("sx", "σ<sub>X</sub> – SD of X", 0.5, 1.5, 1, 0.1, COLS[["X"]]),
slider("sy", "σ<sub>Y</sub> – SD of Y", 0.5, 1.5, 1, 0.1, COLS[["Y"]]),
slider("mx", "μ<sub>X</sub> – mean of X", -1, 1, 0, 0.1, COLS[["X"]]),
slider("my", "μ<sub>Y</sub> – mean of Y", -1, 1, 0, 0.1, COLS[["Y"]]),
tags$hr(),
h5("Condition on X"),
slider("at", "Slice at X =", -2.5, 2.5, 0, 0.1, COND_COL,
animate = animationOptions(interval = 400, loop = TRUE))
),
div(class = "main",
div(class = "card c3d",
div(class = "hdr",
h4("Joint density f(x, y)"),
div(class = "opts",
checkboxInput("show_slice", "Show slice", TRUE),
checkboxInput("show_cef", "E[Y | X = x] line", FALSE),
tags$a(id = "reset3d", "reset view"),
span("drag to rotate · scroll to zoom")
)
),
div(id = "bvn3d")
),
div(class = "card c2d",
div(class = "hdr", uiOutput("cond_title", inline = TRUE),
span(class = "model", "the orange slice, rescaled to area 1")),
uiOutput("stats"),
div(class = "plot2d", plotOutput("cond_plot", height = "100%"))
)
)
),
tags$script(HTML(js))
)
draw_cond <- function(s) {
u <- seq(-LIM_2D, LIM_2D, length.out = 500)
fm <- dnorm(u, s$m_marg, s$s_marg)
fc <- dnorm(u, s$m_cond, s$s_cond)
mc <- COLS[["Y"]]
par(mar = c(2.7, 3.4, 0.5, 0.6), mgp = c(1.6, 0.5, 0), tcl = -0.25, las = 1,
col.axis = "#555", family = "sans")
plot(NA, xlim = c(-LIM_2D, LIM_2D), ylim = c(0, 1.15 * max(fm, fc)), xaxs = "i", yaxs = "i",
xlab = "y", ylab = "density", bty = "l", fg = "#999")
polygon(c(u, rev(u)), c(fm, 0 * u), col = adjustcolor(mc, 0.15), border = NA)
polygon(c(u, rev(u)), c(fc, 0 * u), col = adjustcolor(COND_COL, 0.25), border = NA)
segments(s$m_marg, 0, s$m_marg, max(fm), col = mc, lty = 3)
segments(s$m_cond, 0, s$m_cond, max(fc), col = COND_COL, lty = 2)
lines(u, fm, col = mc, lwd = 2.5)
lines(u, fc, col = COND_COL, lwd = 3)
}
server <- function(input, output, session) {
pars <- reactive(list(mx = input$mx, my = input$my, sx = input$sx, sy = input$sy,
rho = input$rho))
sl <- reactive(slice_of(pars(), input$at))
observe({
p <- pars()
s <- sl()
z <- outer(G, G, function(y, x) dbvn(x, y, p)) # rows are y, as plotly wants
peak <- 1 / (2 * pi * p$sx * p$sy * sqrt(1 - p$rho^2))
# Conditional-mean line E[Y | X = x] across the floor, kept inside the box
cef <- p$my + s$slope * (G - p$mx)
keep <- abs(cef) <= LIM
session$sendCustomMessage("bvn", list(
grid = G,
z = signif(z, 4),
peak = peak,
margX = dnorm(G, p$mx, p$sx),
margY = dnorm(G, p$my, p$sy),
# The box fits the surface and the marginals, with a floor so it holds
# still over the usual range of the sliders
zmax = max(0.45, 1.08 * max(peak, dnorm(0) / min(p$sx, p$sy))),
slice = list(show = input$show_slice, at = s$at, z = dbvn(s$at, G, p), mean = s$m_cond),
cef = list(show = input$show_cef, x = G[keep], y = cef[keep])
))
})
# The title doubles as the plot's legend: conditional in orange, marginal in
# Y's blue
output$cond_title <- renderUI({
h4(span(class = "c", sprintf("Y given X = %.1f", sl()$at)), " against ",
span(class = "m", "Y on its own"))
})
output$stats <- renderUI({
s <- sl()
r <- s$s_cond / s$s_marg
div(class = "fit",
span(class = "m", "SD(Y) ", tags$b(sprintf("%.2f", s$s_marg))),
span(class = "c", sprintf("SD(Y | X = %.1f) ", s$at),
tags$b(sprintf("%.2f", s$s_cond)),
span(class = "meter",
div(style = sprintf("width:%.1f%%;background:%s", 100 * r, COND_COL)))),
span("ratio √(1 − ρ²) ", tags$b(sprintf("%.2f", r))),
span(class = "c", sprintf("E[Y | X = %.1f] ", s$at),
tags$b(sprintf("%.2f", s$m_cond)))
)
})
output$cond_plot <- renderPlot(draw_cond(sl()), res = 96)
}
shinyApp(ui, server)