推文內(nèi)容來自于鏈接
https://www.andrewheiss.com/blog/2021/12/18/bayesian-propensity-scores-weights/
這個博文里的內(nèi)容還挺多的逾雄,我們只關(guān)注其中關(guān)于頻率分布直方圖的實現(xiàn)代碼泽裳。
讀取數(shù)據(jù)集
nets_with_weights<-read.csv("nets_with_weights.csv")
準(zhǔn)備作圖配色
isfahan <- MetBrewer::met.brewer("Isfahan1")
length(isfahan)
isfahan[1]
這里用到的配色包是 https://github.com/BlakeRMills/MetBrewer 這個用到的都是博物館里的油畫的配色寇漫,挺有意思的嚼隘,大家可以試試
使用ggplot2作圖
這里頻率分布直方圖用到的是geom_histogram()
函數(shù),這里的代碼多了一個weight
參數(shù)澈蝙,暫時沒有想明白這個參數(shù)起到什么作用
還遇到一個新函數(shù)colorspace::lighten()
操作顏色祝闻,看幫助文檔是是顏色更亮仪缸。做一個散點圖試試效果
library(ggplot2)
library(patchwork)
p1<-ggplot()+
geom_point(aes(x=1,y=1),size=50,color="darkgreen")
p2<-ggplot()+
geom_point(aes(x=1,y=1),size=50,
color=colorspace::lighten("darkgreen",0.9))
p1+p2
頻率分布直方圖
ggplot() +
geom_histogram(data = filter(nets_with_weights, net_num == 1),
bins = 50,
aes(x = propensity, weight = iptw),
fill = colorspace::lighten(isfahan[2], 0.35),
color="white")
如果要倒過來加一個負(fù)號就可以了
ggplot() +
geom_histogram(data = filter(nets_with_weights, net_num == 1),
bins = 50,
aes(x = propensity, weight = iptw),
fill = colorspace::lighten(isfahan[2], 0.35),
color="white")+
geom_histogram(data = filter(nets_with_weights, net_num == 0),
bins = 50, aes(x = propensity, weight = iptw,
y = -..count..),
fill = colorspace::lighten(isfahan[6], 0.35),
color="white")
添加文本注釋
ggplot() +
geom_histogram(data = filter(nets_with_weights, net_num == 1),
bins = 50,
aes(x = propensity, weight = iptw),
fill = colorspace::lighten(isfahan[2], 0.35),
color="white")+
geom_histogram(data = filter(nets_with_weights, net_num == 0),
bins = 50, aes(x = propensity, weight = iptw,
y = -..count..),
fill = colorspace::lighten(isfahan[6], 0.35),
color="white")+
geom_histogram(data = filter(nets_with_weights, net_num == 1),
bins = 50, aes(x = propensity),
fill = isfahan[2],color="white") +
geom_histogram(data = filter(nets_with_weights, net_num == 0),
bins = 50, aes(x = propensity, y = -..count..),
fill = isfahan[6],
color="white")+
annotate(geom = "label",
x = 0.8, y = 70,
label = "Treated (actual)",
fill = isfahan[2],
color = "white", hjust = 1) +
annotate(geom = "label", x = 0.8,
y = 90, label = "Treated (IPTW pseudo-population)",
fill = colorspace::lighten(isfahan[2], 0.35),
color = "white", hjust = 1) +
annotate(geom = "label", x = 0.8, y = -60,
label = "Untreated (actual)",
fill = isfahan[6],
color = "white", hjust = 1) +
annotate(geom = "label",
x = 0.8, y = -80,
label = "Untreated (IPTW pseudo-population)",
fill = colorspace::lighten(isfahan[6], 0.35),
color = "white", hjust = 1)
對細(xì)節(jié)的一些調(diào)整
ggplot() +
geom_histogram(data = filter(nets_with_weights, net_num == 1),
bins = 50,
aes(x = propensity, weight = iptw),
fill = colorspace::lighten(isfahan[2], 0.35),
color="white")+
geom_histogram(data = filter(nets_with_weights, net_num == 0),
bins = 50, aes(x = propensity, weight = iptw,
y = -..count..),
fill = colorspace::lighten(isfahan[6], 0.35),
color="white")+
geom_histogram(data = filter(nets_with_weights, net_num == 1),
bins = 50, aes(x = propensity),
fill = isfahan[2],color="white") +
geom_histogram(data = filter(nets_with_weights, net_num == 0),
bins = 50, aes(x = propensity, y = -..count..),
fill = isfahan[6],
color="white")+
annotate(geom = "label",
x = 0.8, y = 70,
label = "Treated (actual)",
fill = isfahan[2],
color = "white", hjust = 1) +
annotate(geom = "label", x = 0.8,
y = 90, label = "Treated (IPTW pseudo-population)",
fill = colorspace::lighten(isfahan[2], 0.35),
color = "white", hjust = 1) +
annotate(geom = "label", x = 0.8, y = -60,
label = "Untreated (actual)",
fill = isfahan[6],
color = "white", hjust = 1) +
annotate(geom = "label",
x = 0.8, y = -80,
label = "Untreated (IPTW pseudo-population)",
fill = colorspace::lighten(isfahan[6], 0.35),
color = "white", hjust = 1) +
geom_hline(yintercept = 0, color = "white", size = 0.25) +
scale_y_continuous(label = abs) +
coord_cartesian(xlim = c(0.1, 0.8), ylim = c(-80, 100)) +
labs(x = "Propensity", y = "Count")+
theme_minimal() +
theme(panel.grid.minor = element_blank(),
plot.background = element_rect(fill = "white", color = NA),
plot.title = element_text(face = "bold"),
axis.title = element_text(face = "bold"),
strip.text = element_text(face = "bold", size = rel(0.8), hjust = 0),
strip.background = element_rect(fill = "grey80", color = NA),
legend.title = element_text(face = "bold"))
歡迎大家關(guān)注我的公眾號
小明的數(shù)據(jù)分析筆記本
小明的數(shù)據(jù)分析筆記本 公眾號 主要分享:1、R語言和python做數(shù)據(jù)分析和數(shù)據(jù)可視化的簡單小例子痴荐;2血柳、園藝植物相關(guān)轉(zhuǎn)錄組學(xué)、基因組學(xué)生兆、群體遺傳學(xué)文獻(xiàn)閱讀筆記难捌;3、生物信息學(xué)入門學(xué)習(xí)資料及自己的學(xué)習(xí)筆記鸦难!