バレンタインなのでStanで絵を描いた

こちら、LOVE Equationです。理系のカップルが送ったりするそうです。

x^2 + \left(y - \sqrt[3]{x^2}\right)^2 = 1

僕も、Stanで描いてみました。

heart.stan

data {
  real<lower=0> lambda;
}

parameters {
  real x;
  real y;
}

model {
  target += -lambda*(x^2.0 + (y - (x^2.0)^(1.0/3.0))^2.0 - 1)^2;
}

heart.R

library(rstan)
library(ggplot2)

pars <- c("x", "y")
warmup <- 1000
iter <- 3000
seed <- 87
chains <- 4
src <- "heart.stan"
lambda <- 100

basefit <- stan(file = src, 
                data = list(lambda = lambda), pars = pars,
                seed = seed, chains = chains,
                warmup = warmup, iter = iter, verbose = FALSE)


d <- NULL
for (lambda in c(1, 10, 100, 1000))
{
  fit <- stan(fit = basefit,
              data = list(lambda = lambda), pars = pars,
              seed = seed, chains = chains,
              warmup = warmup, iter = iter, verbose = FALSE)
  tmp <- as.data.frame(fit)
  tmp$lambda <- lambda
  d <- rbind(d, tmp)
}


ggplot(d, aes(x, y)) + 
  geom_point(color = "red") + 
  facet_wrap(~lambda) + coord_fixed() + theme_bw()