ggplot vs Mathematica

หลังๆมานี้ผมใช้งาน ggplot เยอะพอสมควรด้วยความที่มัน plot แล้วลงสีค่าตัวแปรที่ผมต้องการได้ง่าย อย่างเช่นในตัวอย่างนี้ จากชุดข้อมูล iris ผม plot จุดโดยแกน x คือ sepal length และ y คือ sepal width แล้วให้ลงสีตาม species เพื่อที่เราจะได้เห็นว่าแต่ล่ะ species แต่กต่างกันอย่างไร ใน ggplot ผมเพียงแค่บอกว่า แกน x y คือค่าอะไร และลงสีโดยอาศัยค่าจากตัวแปลใด ซึ่งมันดูง่ายมากครับ

กลับมาที่ตัว Mathematica ที่ผมใช้อยู่ประจำต้องบอกเลยครับการจะทำแบบที่ ggplot ทำได้นั้นมันก็ไม่ยากเลยแต่ต้องเขียน code เพิ่มนิดหน่อยในการที่ให้ได้ list ของจุดในแต่ล่ะ species แล้ว plot รวมกันหรือจะแยก plot คนล่ะสีแล้วใช้คำสั่ง Show เอามาซ้อนๆกันได้ ซึ่งผมเองดูแล้วมันก็ค่อนข้างยุ่งพอสมควร แต่ผมก็มาพบด้วยความบังเอิญว่า เราสามารถใช้ Dataset กับคำสั่งพวก plot กราฟ ต่างๆได้ และมันก็ลงสีให้อัตโนมัติด้วย เช่นเพื่อให้ได้กราฟแบบที่ ggplot ทำได้ผมก็แค่ group มันด้วย species ก่อนแล้วเลือก แกนที่จะนำมา plot เช่น

มันทำให้ผมประหลาดใจมาครับที่ Mathematica ทำแบบนี้ได้ตั้งนานแล้ว และผมเองก็เพิ่งรู้ 555

มิกกี้เม้าส์ใน R

วาด Mickey mouse ใน R

library(showtext)
library(ggplot2)

font_add("TH Mali Grade 6", 
"PATH_to_/TH Mali Grade6.ttf")

showtext_auto()
# Define the functions for x(t) and y(t)
x_func <- function(t) {
  -119/9 * sin(1/10 - 25 * t) - 945/43 * sin(5/11 - 23 * t) - 501/8 * sin(1/47 - 10 * t) - 
    801/20 * sin(13/10 - 8 * t) - 607/10 * sin(7/10 - 7 * t) - 1975/9 * sin(3/7 - 2 * t) + 
    4655/8 * sin(t + 29/10) + 773/7 * sin(3 * t + 76/17) + 1138/9 * sin(4 * t + 1/3) + 
    619/4 * sin(5 * t + 5/8) + 3355/26 * sin(6 * t + 13/10) + 595/9 * sin(9 * t + 28/11) + 
    905/9 * sin(11 * t + 40/11) + 149/10 * sin(12 * t + 7/11) + 153/5 * sin(13 * t + 2) + 
    127/8 * sin(14 * t + 13/6) + 225/8 * sin(15 * t + 3/2) + 49/4 * sin(16 * t + 25/7) + 
    133/9 * sin(17 * t + 15/8) + 301/15 * sin(18 * t + 31/9) + 232/7 * sin(19 * t + 11/7) + 
    438/19 * sin(20 * t + 10/7) + 67/3 * sin(21 * t + 19/10) + 114/11 * sin(22 * t + 38/37) + 
    153/11 * sin(24 * t + 10/9) + 153/8 * sin(26 * t + 11/12) + 443/34 * sin(27 * t + 23/14) + 
    39/7 * sin(28 * t + 45/11) + 32/3 * sin(29 * t + 4/11) + 110/13 * sin(30 * t + 3/8) + 
    41/7 * sin(31 * t + 70/23) + 39/7 * sin(32 * t + 13/9)
}

y_func <- function(t) {
  -86/19 * sin(6/7 - 30 * t) - 89/8 * sin(2/3 - 29 * t) - 687/49 * sin(4/9 - 24 * t) - 
    294/11 * sin(11/8 - 20 * t) - 135/4 * sin(10/11 - 17 * t) - 71/7 * sin(1/6 - 15 * t) - 
    226/5 * sin(3/8 - 14 * t) - 295/8 * sin(1/8 - 11 * t) - 220/7 * sin(17/16 - 10 * t) - 
    708/7 * sin(1/6 - 9 * t) - 2017/8 * sin(5/4 - 2 * t) - 1240/3 * sin(11/10 - t) + 
    104/15 * sin(18 * t) + 2216/7 * sin(3 * t + 39/11) + 601/3 * sin(4 * t + 18/11) + 
    468/11 * sin(5 * t + 34/11) + 127/7 * sin(6 * t + 27/7) + 457/8 * sin(7 * t + 25/6) + 
    201/4 * sin(8 * t + 23/7) + 71/3 * sin(12 * t + 4/7) + 340/7 * sin(13 * t + 39/11) + 
    216/7 * sin(16 * t + 27/13) + 31/5 * sin(19 * t + 24/25) + 8/5 * sin(21 * t + 1/7) + 
    217/11 * sin(22 * t + 23/14) + 34/3 * sin(23 * t + 19/8) + 157/9 * sin(25 * t + 3/10) + 
    69/4 * sin(26 * t + 2/3) + 29/5 * sin(27 * t + 15/7) + 58/7 * sin(28 * t + 29/8) + 
    95/11 * sin(31 * t + 19/8) + 34/9 * sin(32 * t + 40/9)
}

# Generate values
t_values <- seq(0, 2 * pi, length.out = 1000)
x_values <- sapply(t_values, x_func)
y_values <- sapply(t_values, y_func)

df <- data.frame(t = t_values, x = x_values, y = y_values)

# Plot using ggplot2
p <- ggplot(df, aes(x = x, y = y)) + 
  geom_path(color = "black", size = 1) + 
  theme_minimal() + 
  geom_text(aes(x = 500, y = 800, label = "กลุ่มผู้ใช้ R"), 
            family = "TH Mali Grade 6", color = "black", size = 26)

print(p)

เรียน R Shiny

สำหรับผู้ที่สนใจจะพัฒนาเว็บและแอปใน R ด้วย Shiny ครับ เนื้อหาที่ผมสอนใน workshop ซึ่งเป็นภาษาไทย ได้ถูกเผยแพร่ไว้ที่นี่ Shiny เบื้องต้น (ammnet-thailand.github.io) ครับ

Merry Christmas and a Happy New Year!!

library(rgl)
open3d()
bg3d("black")
t <- seq(0,8*pi, by = 0.01);
x <- c(t*cos(t),t*cos(t+0.5*pi),t*cos(t+pi),t*cos(t+1.5*pi))
y <- c(t*sin(t),t*sin(t+0.5*pi),t*sin(t+pi),t*sin(t+1.5*pi))
z <- -c(t,t,t)
plot3d(x, y, z, col = rainbow(1000),box=FALSE,axes=FALSE,xlab="",ylab="",zlab="")
text3d(1,8, text = "Merry Christmas and", adj = 0.5, color = "green")
text3d(1,4, text = "a Happy New Year!", adj = 0.5, color = "green")
points3d(runif(100,min=-25,max = 25),runif(100,min=-25,max = 25),-runif(100,max = 25),color="white")

Kung package แพ็คเกจกุ้ง

แพ็คเกจกุ้งเป็นแพ็คเกจภาษา R ที่ผมเขียนขึ้นมาสำหรับใช้สร้าง Shiny App สำหรับโค้ดของโมเดลที่สร้างจากสมการอนุพันธ์หรือ ODE ที่ใช้ตัว solver จากแพ็คเกจที่ชื่อ deSolve ครับ

ไอเดียก็มีเพียงว่าจากสมการ ode ที่สร้างขึ้นมาด้วยถ้าตัวแปรในสมการมีการเปลี่ยนแปลงค่า ผลลัพท์ที่ออกมามันจะเปลี่ยนแปลงไปอย่างไร การเขียนใน R นั้นเราก็อาจจะใช้คำสั่ง manipulate มาช่วยได้แต่มันก็ยังมีข้อจำกัดเรื่องจำนวนของตัวแปรที่สามารถเรียกมาใช้ได้ เพราะมันขึ้นอยู่กับขนาดหน้าจอของตัว plot ผมเลยคิดว่าถ้าใช้ shiny มันน่าจะดูดีกว่าและสามารถใส่ตัวแปรได้มากกว่า ขึ้นกับการออกแบบ ui แต่การเขียน shiny ก็ไม่ใช่ว่าจะง่ายสำหรับคนเริ่มเรียนรู้ R ผมเลยคิดว่าน่าจะเขียนแพ็คเกจที่สามารถช่วยเรื่องนี้ได้

แพ็คเกจที่เขียนก็พยายามทำให้มันใช้ง่ายมากที่สุดจาก code ที่มีอยู่แล้ว โดยผู้ใช้เพิ่มคำสั่งไม่มากก็สามารถที่จะสร้าง shiny application ได้แล้ว ตัวอย่างเช่น ถ้าผู้ใช้มี code อยู่แล้วตามนี้ ซึ่งมันจะคำนวณและ plot กราฟให้ตามค่าของตัวแปรพารามิเตอร์ที่กำหนดไว้ เช่น gamma = 0.14286 และ beta = 0.6

library(deSolve)

init <- c(S = 1-1e-6, I = 1e-6, R= 0.0)
times <- seq(0, 70, by = 1)

sir <-function(time, state, parameters) {
  with(as.list(c(state, parameters)), {
    dS <- -beta*S*I
    dI <- beta*S*I - gamma*I
    dR <- gamma*I
    
    return(list(c(dS, dI, dR)))
  })
}

parameters <- c(
  gamma = 0.14286, 
  beta = 0.6
)

out <- as.data.frame(ode(y = init, times = times, func = sir, parms = parameters))

matplot(times, out[,c("S","I","R")], type = "l", xlab = "Time", ylab = "Susceptibles and Recovereds", main = "SIR Model", lwd = 1, lty = 1, bty = "l", col = 2:4)
legend(40, 0.7, c("Susceptibles", "Infecteds", "Recovereds"), pch = 1, col = 2:4)

ถ้าผู้ใช้ต้องการดูว่ากราฟจะเปลี่ยนแปลงอย่างไรถ้าตัวแปร gamma หรือ beta มีการเปลี่ยน ก็เพียงเพิ่ม code ลงไปในส่วนต่างๆ ตามนี้

!Start

library(deSolve)

init <- c(S = 1-1e-6, I = 1e-6, R= 0.0)
times <- seq(0, 70, by = 1)

sir <-function(time, state, parameters) {
  with(as.list(c(state, parameters)), {
    dS <- -beta*S*I
    dI <- beta*S*I - gamma*I
    dR <- gamma*I

    return(list(c(dS, dI, dR)))
  })
}

!Parameters
parameters <- c(
  gamma = 0.14286,
  beta = 0.6
)

!ODECMD
out <- as.data.frame(ode(y = init, times = times, func = sir, parms = parameters))

!PostProcess

!Plots
matplot(times, out[,c("S","I","R")], type = "l", xlab = "Time", ylab = "Susceptibles and Recovereds", main = "SIR Model", lwd = 1, lty = 1, bty = "l", col = 2:4)
legend(40, 0.7, c("Susceptibles", "Infecteds", "Recovereds"), pch = 1, col = 2:4)

!Controls
sliderInput("beta","beta", min = 0,max = 10,step = 0.001,value = 1.4),
sliderInput("gamma","gamma", min = 0,max = 1,step = 0.001,value = 0.14)

!End

จะเห็นว่ามี keywords ที่เพิ่มเข้าไปดังนี้

!Start – เป็นตัวบอกเริ่มต้น เราสามารถบอกรายละเอียดของตัวแปรอย่าง เวลา และค่าเริ่มต้นของ state หรือ compartment ต่างๆและ function ของระบบสมการหลักที่จะแก้ด้วย deSolve::ode ใน keywords นี้ได้เลยครับ

!Parameters – เป็นตัวแปร parameters ของโมเดล ode ครับ ซึ่งต้องกำหนดเป็น vector ตามตัวอย่างนี้เลยครับ

!ODECMD – สำหรับคำสั่งที่ต้องใช้คำสั่ง ode จาก deSolve ครับ สำหรับตอนนี้ต้องให้ตัวแปร output ชื่อ out นะครับ

!PostProcess – ถ้ามีการคำนวณค่าอะไรบางอย่างจาก out ก็ให้ใส่ใน keyword นี้ได้เลยครับ

!Controls – สำหรับใส่ sliders ของตัวแปรที่สนใจครับ

!End – สำหรับบอกว่าจบแล้ว 🙂

โดยหลังจากที่ใส่ keywords และsave เป็นไฟล์ใหม่แล้ว เช่นชื่อ mysystem.R เราสามารถสร้าง shiny app ได้ด้วยการรันคำสั่งนี้ครับ

runSystem(‘mysystem.R’)

ผลที่ได้ก็จะประมาณนี้ครับ

ถ้าสนใจอยากใช้งานก็ทำการติดตั้งได้โดยพิมพ์

devtools::install_github('slphyx/Kung')

ดูเพิ่มเติมที่ https://github.com/slphyx/Kung

%d bloggers like this: