520怎么办？用R画心形函数｜来点理工式浪漫！！

1. 前言

不务正业一下，笛卡尔坐标系应该大家都知道，心形函数应该也略有听闻，

摘抄一下GPT：

关于心形函数的爱情故事虽然更多是浪漫的传说和文学创作，而非确凿的历史事实，但其中一个广为流传且感人的故事涉及到17世纪法国数学家笛卡尔（René Descartes）。

传说中的爱情故事

据传说，笛卡尔爱上了一位美丽的公主（有时被称为“伊丽莎白公主”）。笛卡尔希望通过一种独特而永恒的方式来表达他对公主的爱。于是，他利用他在数学上的天才，设计了一个数学方程，该方程描绘出一个完美的心形曲线。这条曲线不仅象征着他对公主的爱，还展示了数学之美。

故事进一步描述，笛卡尔将这个心形曲线方程写在纸上，并送给公主作为表达他爱意的礼物。公主被笛卡尔的聪明才智和浪漫情怀所打动，两人最终坠入爱河。

尽管这个故事没有历史记录支持，但它在许多浪漫主义者和数学爱好者中广为流传，成为一个美丽的象征。这个故事不仅仅展示了笛卡尔的数学才华，还强调了数学在表达情感和美感方面的独特作用。

2.公式

2.1 极坐标图

现在我们进入数学表示。在笛卡尔坐标系中，点表示为 (x,y)，其中 x 轴和 y 轴彼此垂直。在极坐标系中，点表示为 (r, θ)，其中 r 是距原点（“半径”）的距离，θ 是与 x 轴正值部分的逆时针角度在笛卡尔图上。

r=1-cosθ

2.2 二维隐式曲线

隐式方程是其中一个变量不作为另一变量的函数给出的方程，通常是因为显式表达一个变量非常困难，或者因为隐式关系更简单；相反，这两个变量都用于定义函数。上面的心形曲线可以在笛卡尔空间中隐式定义为

(x^2+y^2-1)^3-x^2y^3=0

2.3 参数曲线

参数方程通过将每个分量x和y（有时还有z）作为t的函数来定义曲线。非常困难甚至不可能表达为显式或隐式方程的曲线通常可以非常简单地表达为参数方程。上面的心形曲线可以表示为参数方程

(x=16\sin^3t, y=13\cos t-5\cos(2t)-2\cos(3t)-cos(4t))

3.代码

先上个青春版的心形函数：

library(ggplot2)
dat <- data.frame(t = seq(0, 2*pi, by = 0.01))
x <-  function(t) 16 * sin(t)^3
y <- function(t) 13*cos(t) - 5*cos(2*t) - 2*cos(3*t) - cos(4*t)

dat$y <- y(dat$t)
dat$x <- x(dat$t)

heart <- ggplot(dat, aes(x,y)) +
  geom_polygon(fill = "red", col = "firebrick", alpha = 0.9) +
  theme_classic()+
  title("r=a(1-sinθ)")

print(heart)

进阶版：首先定义了一个心形曲线的数学公式，并计算了曲线上的点，然后生成一个数据框以随机填充心形，并赋予每个点不同的颜色和大小。接着，使用ggplot2包绘制并保存心形图为PNG图像。最后，通过循环逐步增加点的数量和位置，生成一个逐渐填充心形并最终漂浮起来的GIF动画，使用saveGIF函数将动画保存为文件。

生成的输入文件格式：

# requires ImageMagick (http://www.imagemagick.org/)
# install.packages("animation", repos = "http://rforge.net", type = "source")
# install.packages("dplyr")
# install.packages("ggplot2")
# install.packages("pryr")

require(animation)
require(dplyr)
require(ggplot2)
require(pryr)
# install.packages("pryr")
#------------------------------------------------------------------------------#
# heart curve formula
heart <- quote((x^2 + y^2 - 1)^3 - x^2 * y^3)

# formula for heart curve at a given x
heart_at_x <- function(x) {
  function(y) eval(substitute_q(heart, list(x = x)), list(y = y))
}

# trace the heart curve
# by evaluating the heart curve formula at each x, then finding the roots of the
# resulting formula in y; e.g. a x==0, find the roots of (y^2 - 1)^3 = 0
# broken up into upper and lower parts (h_y1 and h_y2)
heart_x <- seq(-1.136, 1.136, 0.001)
heart_y_lower <- sapply(heart_x, function(x) uniroot(heart_at_x(x), c(-2, 0.6))$root)
heart_y_upper <- sapply(heart_x, function(x) uniroot(heart_at_x(x), c(0.6, 2))$root)

# put together data frame
heart_df <- data.frame(x = rep(heart_x, 2), 
                       y = c(heart_y_lower, heart_y_upper))

# show outline
with(heart_df, plot(x, y))

# create a data frame with one row per x, so we can fill in the heart
heart_df_minmax <- data.frame(x = heart_x,  
                              y_min = heart_y_lower, 
                              y_max = heart_y_upper)

set.seed(20150214)

# fill in the heart by generating random deviates at each x 
# and rejecting those that fall outside the heart curve
heart_full <- apply(heart_df_minmax, 
                    1, 
                    function(w) {
                      x <- w["x"]
                      y_min = w["y_min"]
                      y_max = w["y_max"]
                      y <- rnorm(2, mean = 0.33)
                      y <- y[between(y, y_min, y_max)]
                      x <- x[any(is.finite(y))]
                      data.frame(x, y, row.names = NULL)
                    })

# change from list to data frame
heart_full <- bind_rows(heart_full)

# add random numbers for color and size
heart_full <- heart_full %>% 
  mutate(z1 = runif(n()), 
         z2 = pmin(abs(rnorm(n())), 3), 
         order = runif(n())) %>%
  arrange(order)
# colnames(heart_full) <- c("20240520"," "," "," ","orders")
#------------------------------------------------------------------------------#

# plot the heart
library(ggplot2)
library(extrafont)
p <- ggplot(heart_full, 
            aes(x, y, color = z1, size = z2)) + 
  geom_point(pch = -1 * as.hexmode(9829)) + 
  scale_color_gradient(limits = c(0, 1), low = "red2", high = "pink") + 
  scale_size(limits = c(0, 3), range = c(0.1, 20)) + 
  theme_bw() +
  labs(title = "吃好喝好 快高長大 開心每一天～",  
       x = "Dear Rose",    
       y = "20240520",    
       color = "sth.strength",       
       size = "sth.size") +      
  theme(plot.title = element_text(family = "SimHei", face = "bold", color = "#6495ed", size = 24, hjust = 0.5),
        axis.title.x = element_text(family = "sans", face = "bold", color = "#ff4500", size = 14),
        axis.title.y = element_text(family = "SimHei", face = "bold", color = "#ff4500", size = 14),
        legend.title = element_text(family = "sans", face = "bold", color = "black", size = 11))

print(p)  
png("valentine.png", 800, 600)
# pdf("dingdingdang.pdf",width = 10,height = 7)
print(p) 
dev.off()

# animated plot
saveGIF({
  fill_steps <- 50  # heart fill-in frames
  float_steps <- 25  # heart float-away frames
  
  for (i in seq(fill_steps + float_steps)) {
    # find the number of hearts to fill in on this step
    num_hearts <- min(i, fill_steps) * nrow(heart_full) / fill_steps
    # once the heart is filled in, make the heart float away
    # by shifting each point up some amount
    if (i > fill_steps) {
      j <- i - fill_steps
      j_scale <- uniroot(function(x) (x * float_steps)^2 - 2.5, c(0, 1))$root
      y_change <- (j_scale * j)^2
      heart_full <- mutate(heart_full, y = y + y_change)
    }
    # plot the heart
    p <- ggplot(heart_full[seq(num_hearts), ], 
                aes(x, y, color = z1, size = z2)) + 
      geom_point(pch = -1 * as.hexmode(9829)) + 
      scale_color_gradient(limits = c(0, 1), low = "red3", high = "pink") + 
      scale_size(limits = c(0, 3), range = c(0.1, 20)) + 
      theme_bw() + 
      coord_cartesian(xlim = c(-1.5, 1.5), ylim = c(-1.25, 1.5))
    print(p)
  }
}, 
movie.name = "valentine.gif", 
interval = 0.1, 
nmax = 30, 
ani.width = 600, 
ani.height = 400)

其实还有三维的，莫得时间折腾了，心意到位即可～