odbcパッケージ経由だとvarcharが255文字に切り詰められる問題

ハマるの2回目なのでメモっておく。

PostgreSQLで“test”という名前のデータベースが存在していて、そこに接続している状態を想定する。

con <- odbc::dbConnect(odbc::odbc(), "test")

PostgreSQLのcharacter varying(varchar)は、長さ指定なしで使うと文字列の長さの制限がなくなる(cf. 8.3. 文字型)。

しかし、これをodbcパッケージ経由で取得すると255文字に切り詰められる(https://github.com/r-dbi/odbc/issues/202 で報告されているのと同じ現象?)。

library(dplyr)
query <- "SELECT (repeat('a', 1000))::varchar as str;"
odbc::dbGetQuery(con, glue::glue("CREATE TABLE x AS {query}"))
odbc::dbGetQuery(con, "SELECT str FROM x;") %>% nchar()
## str 
## 255

これは、TEXT型にキャストするととりあえず回避できる。

odbc::dbGetQuery(con, "SELECT CAST(str AS TEXT) FROM x;") %>% nchar()
##  str 
## 1000

dplyrの場合は、as.character()を使う。

tbl(con, "x") %>% mutate(str = as.character(str)) %>% pull() %>% nchar()
## [1] 1000

as.character()はTEXT型へのキャストをするSQLを生成するので、やっていることは同じ。

dbplyr::translate_sql(as.character(x))
## <SQL> CAST("x" AS TEXT)

「マンガでわかる統計学[因子分析編]」の主成分分析をRでやる

やったのでメモっておく。

「第4章 主成分分析」より。

主成分と主成分得点を求める

Step0. データ入力

library(tibble)
ramen <- tribble(
  ~店舗,    ~, ~, ~スープ,
  "二郎",       2, 4, 5,
  "夢田屋",     1, 5, 1,
  "地回",       5, 3, 4,
  "なのはな",   2, 2, 3,
  "花ぶし",     3, 5, 5,
  "昇辰軒",     4, 3, 2,
  "丸蔵らあめん", 4, 4, 3,
  "海楽亭",     1, 2, 1,
  "なるみ家",   3, 3, 2,
  "奏月",       5, 5, 3
)

Step1. 変数ごとに基準化する

library(dplyr)

ramen %>% 
  mutate_if(is.numeric, scale) ->
  ramen_scaled
ramen_scaled
## # A tibble: 10 x 4
##    店舗             麺     具  スープ
##    <chr>         <dbl>  <dbl>   <dbl>
##  1 二郎         -0.671  0.341  1.45  
##  2 夢田屋       -1.34   1.19  -1.31  
##  3 地回          1.34  -0.511  0.759 
##  4 なのはな     -0.671 -1.36   0.0690
##  5 花ぶし        0      1.19   1.45  
##  6 昇辰軒        0.671 -0.511 -0.621 
##  7 丸蔵らあめん  0.671  0.341  0.0690
##  8 海楽亭       -1.34  -1.36  -1.31  
##  9 なるみ家      0     -0.511 -0.621 
## 10 奏月          1.34   1.19   0.0690

Step2. 相関行列を求める

ramen_scaled %>% 
  select(-店舗) %>% 
  cor ->
  ramen_cormat
ramen_cormat
##               麺        具    スープ
## 麺     1.0000000 0.1905002 0.3600411
## 具     0.1905002 1.0000000 0.3004804
## スープ 0.3600411 0.3004804 1.0000000

Step3. 固有値固有ベクトルを求める

ramen_decomp <- eigen(ramen_cormat)
ramen_decomp
## eigen() decomposition
## $values
## [1] 1.5728539 0.8140083 0.6131378
## 
## $vectors
##            [,1]       [,2]       [,3]
## [1,] -0.5715110  0.6044710  0.5549685
## [2,] -0.5221161 -0.7896069  0.3223595
## [3,] -0.6330639  0.1055260 -0.7668731

固有ベクトルは定数倍の任意性があるので、テキストと一緒にならない可能性があり、実際なっていない。

Step4

library(ggplot2)
library(ggrepel)
target <- 1:2
target_vec <- ramen_decomp$vectors[, target]
as.data.frame(target_vec) %>% 
  mutate(item = c("麺", "具", "スープ")) %>% 
  ggplot(aes(x = -V1, y = -V2)) + # テキストに合わせるために反転
  geom_point() +
  geom_text_repel(aes(label = item), family = "IPAexGothic") +
  lims(x = c(0, 1), y = c(-1, 1)) +
  geom_vline(xintercept = 0, col = "gray") +
  geom_hline(yintercept = 0, col = "gray") +
  theme_minimal()

f:id:Rion778:20181003211033p:plain

Step5. 第一主成分と第二主成分を理解する。

Step6 各個体の第一主成分と第二主成分における座標、つまり各個体の第一主成分得点と第二主成分得点を求める。

ramen_scaled %>% 
  rowwise() %>% 
  mutate(
    score1 = -sum(c(,, スープ) * ramen_decomp$vectors[, 1]),
    score2 = -sum(c(,, スープ) * ramen_decomp$vectors[, 2])
  )
## Source: local data frame [10 x 6]
## Groups: <by row>
## 
## # A tibble: 10 x 6
##    店舗             麺     具  スープ score1 score2
##    <chr>         <dbl>  <dbl>   <dbl>  <dbl>  <dbl>
##  1 二郎         -0.671  0.341  1.45    0.712  0.522
##  2 夢田屋       -1.34   1.19  -1.31   -0.974  1.89 
##  3 地回          1.34  -0.511  0.759   0.980 -1.29 
##  4 なのはな     -0.671 -1.36   0.0690 -1.05  -0.678
##  5 花ぶし        0      1.19   1.45    1.54   0.789
##  6 昇辰軒        0.671 -0.511 -0.621  -0.277 -0.744
##  7 丸蔵らあめん  0.671  0.341  0.0690  0.605 -0.144
##  8 海楽亭       -1.34  -1.36  -1.31   -2.31  -0.127
##  9 なるみ家      0     -0.511 -0.621  -0.660 -0.338
## 10 奏月          1.34   1.19   0.0690  1.43   0.124

Step7 第一主成分と第二主成分に基づいてプロットする。

ramen_scaled %>% 
  rowwise() %>% 
  mutate( # テキストに合わせるために-1を掛けて反転させる
    score1 = -sum(c(,, スープ) * ramen_decomp$vectors[, 1]),
    score2 = -sum(c(,, スープ) * ramen_decomp$vectors[, 2])
  ) %>% 
  ggplot(aes(x = score1, y = score2)) +
  geom_point() +
  geom_text_repel(aes(label = 店舗), family = "IPAexGothic") +
  geom_vline(xintercept = 0) +
  geom_hline(yintercept = 0) +
  theme_minimal()

f:id:Rion778:20181003211155p:plain

分析結果の精度を確認する

寄与率を計算。

ramen_decomp$values / 3
## [1] 0.5242846 0.2713361 0.2043793

累積寄与率を計算

cumsum(ramen_decomp$values) / 3 * 100
## [1]  52.42846  79.56207 100.00000

同じことをRでやる

Rで、というかprcompで。

ramen %>%
  select(-店舗) %>% 
  prcomp(scale = TRUE) -> result
result
## Standard deviations (1, .., p=3):
## [1] 1.2541347 0.9022241 0.7830312
## 
## Rotation (n x k) = (3 x 3):
##              PC1        PC2        PC3
## 麺     0.5715110 -0.6044710  0.5549685
## 具     0.5221161  0.7896069  0.3223595
## スープ 0.6330639 -0.1055260 -0.7668731
summary(result)
## Importance of components:
##                           PC1    PC2    PC3
## Standard deviation     1.2541 0.9022 0.7830
## Proportion of Variance 0.5243 0.2713 0.2044
## Cumulative Proportion  0.5243 0.7956 1.0000

先程の計算と一致していることを確認(※主成分に定数倍の差はある)。

標準偏差固有値平方根に一致する。

sqrt(ramen_decomp$values)
## [1] 1.2541347 0.9022241 0.7830312
ramen_decomp
## eigen() decomposition
## $values
## [1] 1.5728539 0.8140083 0.6131378
## 
## $vectors
##            [,1]       [,2]       [,3]
## [1,] -0.5715110  0.6044710  0.5549685
## [2,] -0.5221161 -0.7896069  0.3223595
## [3,] -0.6330639  0.1055260 -0.7668731

結果はbiplotで可視化できる。

biplot(result, family = "IPAexGothic", xlabs = ramen$店舗)

f:id:Rion778:20181003211949p:plain

ところで、ggplotを使いたいと思ってggfortifyggbiplotfactoextraの3つのパッケージを試してみたのだが、現状日本語が問題なく使えるのはggfortifyだけのようだった。

ggfortifyはパッケージロード後にautoplotを通じてプロットをするが、引数を確認するには?ggfortify::ggbiplotを見る必要があった。

library(ggfortify)
autoplot(
  result2, 
  loadings = TRUE, 
  loadings.label = TRUE, 
  loadings.label.family = "IPAexGothic",
  loadings.label.repel = TRUE,
  label = TRUE,
  label.family = "IPAexGothic",
  label.repel = TRUE,
  scale = 0
) +
  theme_bw()

f:id:Rion778:20181003212928p:plain

農薬の名前はポアソン分布に従っている

気がするのでちょっと見て欲しい。

まず、rvestを使ってFAMICから農薬登録情報のcsvをダウンロードして保存する。URLがちょいちょい変わるのでxpathとかで抽出したほうが良い。

library(rvest)
library(dplyr)
url <- "http://www.acis.famic.go.jp/ddata/index2.htm"
read_html(url) %>%
  html_node(xpath = '//*[@id="mainArea"]/ul/li[1]/ul/li/a') %>% 
  html_attr(name = "href") %>% 
  url_absolute(base = url) %>% 
  download.file("kihon.zip")
kihon <- unzip("kihon.zip") %>% 
  readr::read_csv(locale = readr::locale(encoding = "CP932")) %>%
  mutate_if(is.character, stringi::stri_trans_nfkc)

次に農薬の名称を文字数に変換する。比較のために他の文字列情報もいくつか一緒に変換しておく。

kihon %>% 
  select_if(is.character) %>% 
  select(-c(濃度, 用途, 剤型名, 登録年月日, 登録の有効期限)) %>%
  tidyr::gather() %>%
  unique() %>% 
  mutate(n = nchar(value)) -> kihon_mod

次にポアソン分布の確率質量を各変数の平均値から計算しておく。

# 平均値
lambda <- kihon_mod %>% 
  group_by(key) %>% 
  summarise(lambda = mean(n, na.rm = TRUE))
# 確率質量
pois <- data.frame(n = 1:50, key = rep(lambda$key, rep(50, nrow(lambda)))) %>% 
  left_join(lambda, by = "key") %>% 
  mutate(theory = dpois(n, lambda))

そしてプロットする。

library(ggplot2)
library(gghighlight)
kihon_mod  %>% 
  ggplot(aes(x = n)) +
  geom_bar(aes(y = ..prop..)) +
  geom_line(data = pois, 
            aes(x = n, y = theory)) +
  gghighlight(key == "農薬の名称") +
  facet_grid(key ~ ., scales = "free_y") +
  lims(x = c(0, 30)) +
  theme_bw(base_family = "IPAexGothic") 

f:id:Rion778:20180902224525p:plain

農薬の名称はランダムに決まっている…?

ところで最初に雑に「農薬の名前」と言ってしまったが、農薬の名前には「農薬の名称(=商品名)」の他に「種類名」と「化学名」があってややこしいのだが、今回ポアソン分布っぽい分布をしているのは商品名だ。

ちなみに今回プロット対象にした項目の実質的な意味は次の通り。

  • 正式名称...販売しているメーカーの名前。
  • 総使用回数における有効成分...使用回数のカウントの際に用いる有効成分名。なんと実際の有効成分名とイコールではなく厄介(分かりやすい例で言うとこちらでは光学異性体をまとめて1つの扱いにしているものがある)。
  • 農薬の種類...種類名。
  • 農薬の名称...商品名。
  • 有効成分...有効成分名。

なぜもっと単純にならないのか…。

歯を抜いたので状況を記録しておく

親知らずを抜いた。3週間ぶり2本目。次はないので、経過を記録しておきたい。随時追記していく所存。

抜歯にいたる経緯

昨年11月頃、歯が痛くなってしまったので歯医者に行った。昔詰め物をしたところがあったのだが、詰め物が取れてしまい、そのままほっといたら虫歯になってしまったのだ。虫歯の治療は割と痛かったし、セラミックの詰め物(というか被せ物)にしたら財布へのダメージも甚大だったし、歯医者には定期的に行くべきだ。

それはともかく、治療にあたってレントゲンを撮影した。すると、下の両奥歯の奥、明らかに角度を間違えた歯が2本ある。親知らずだ。横向きに生えている親知らずは抜くのが大変だとか聞いたことがあるが、抜かないといけないのだろうか。

「んーまあ左はちょっと頭出てるし、右もほぼ出てるようなもんだし、抜けるでしょう」

抜くことは確定らしい。

横向きに圧力がかかって歯並びに影響しているし、ちょっと出てるから虫歯のリスクも高く、さっさと抜いておくべきとのこと。

何か軽い感じで言ってるし、大した事ないのではないか、その時はそう思っていた。

それからしばらくして虫歯の治療が完了し、親知らずの抜歯の話になった。

「あっでは大学病院に紹介状書きますのでー」

そうだ。横向きの親知らずがそんなに簡単に抜けるわけがない。私は大学病院送りになった。

1本目の抜歯

予約した日に大学病院へ行くと、学生らしき若者に体温と血圧と脈拍を取られる。大学病院っぽさだ。誰だって最初は初心者なのだから仕方がないのだが、一向に脈が取れず、すごく不安感を煽られる。体温は37.8℃あった。どうも極度に緊張したり不安になったりすると発熱する体質らしいと気付いたのはこの時だった。

その日はレントゲン撮影やらなんやらを行った。レントゲン装置が顔の周りをぐるっと回る間、何故か装置からエリーゼのためにが流れる。誰だ無駄な機能を付けたのは。不意打ちはやめろ。そしてその日は日程調整をしただけで抜歯はしなかった。壮絶な決意で行ったのに。両方一度に抜くのではなく、一本抜いて落ち着いたらもう一本抜く、という方針で行くらしい。

日を改めて1本目の抜歯を行った。その日測定してもらった体温は37.6℃だった。親知らずの抜歯はちょっとした手術になるので、リスクの説明とか同意書へのサインとか諸々の儀式があり、不安感を煽られる。

「抜いてる時は最初麻酔でチクッとするだけで痛くないですが、その後は確実に痛いです。3日目くらいが腫れのピークです。」

救いはないのか。

そして手術が始まる。若干想定外だったのは「最初麻酔でチクッと」が1回ではなく3回くらい点だ。割とブスブス刺された気がする。さほど痛くはなかったが。麻酔が効いてしまうと確かに痛みはない。顔には覆いがされているので、後は音で想像するしかない。「水が出ますよ―」と何度か言われたが、正直水が出ているのかどうか、何のために水が出ているのかも良くわからなかった(多分削りカス的なものをアレしたりするためだろうけど想像するのはやめておいた)。

どれくらいの時間がたっただろうか。後で確認したら30分くらいだった。

「抜けましたからねー」

抜けたらしい。左の親知らずだったものは3分割されていた。

「持って帰る?」

いらない。

抜歯後の経過

痛み止めに処方してもらったロキソニンが効いてる間はさほど痛くはない(無痛ではないが、強い不快感があるほどでもない)が、ロキソニンが切れるとすぐ痛くなるという状態が2〜3日は続いた。4日目くらいに血餅が取れてビビったが、これくらい時間が経過していれば大丈夫らしい。

1週間後に再度大学病院に行って抜糸した。抜歯と違って秒単位で終わったが、意外と痛かった。もう一度あると思うと若干不安。この時2本目の抜歯の予約をした。(当然抜くよな)という圧を感じた。選択の余地は無かった。

2週間くらいすると違和感も大分なくなった。ただ穴はあいているので食べ物が入ったりして若干不快感はある。この穴は骨が再生して埋まるらしいが、数ヶ月くらいかかるらしい。

この間、割と大変だったような気もするのだが、喉元過ぎれば何とやらで、細かい部分は忘れてしまった。次回はちゃんと記録したいと思った。


2本目の抜歯

そして2本目の抜歯だ。前回の反省を活かしてこまめに経過を記録していきたい。

抜歯前日まで

前回抜歯時に調査した結果として、抜歯後の早い段階で血餅(かさぶた)が取れてしまうとドライソケットと呼ばれる状態になり、激痛に襲われるという情報を得ていた。しかも下の親知らずはドライソケットのリスクが高いらしい。そこで、少しでも可能性を減らすべく、食糧としてカロリーメイトのゼリータイプを買い込んだ。昔喉をやられたときに1週間くらいこれだけで生きたことがある。個人的にはりんご味が好きだ。

抜歯当日(2018/05/11)

2本目だろうが何本目だろうが嫌なものは嫌なので、めっちゃ嫌だなと思いながら大学病院へ向かう。

体温と血圧と脈拍を測定してもらう。やはり学生、スムーズには行かない。こういう実習も大切と思いつつも不安感が高まる。体温は37.3℃だった。

麻酔をして処置が始まる。前回30分くらいだったので、1800まで数えれば終るはずだと思い、数を数えることに集中する。やたら「破骨鉗子」という単語が聞こえる。何か破片を引っ張り出してるのだろうかと思うが想像するのはやめておく。

「抜けましたからねー」

と言われたのは320くらいまで数えた時だった。きちんと30分経過していたし、30分くらい経ってるよなという感じもあった。人は緊張するとまともに数字を数えられなくなるのだという知見を得た。

右の親知らずだったものは8つくらいの破片に分割されていた。なにか掴んで引っ張り出してるっぽかったのはこれか。

「持って帰る?」

いらない。

ロキソニンは家に在庫があったので、今回は化膿止めだけが処方された。あれ、前回はうがい薬的なものが出た気が…。

不安感を抱えつつその日はそのまま出社し、そっと仕事をした(家より会社のが近いのだ)。

やはりロキソニンが効いてるうちは何ということは無いのだが、薬が切れてくると痛い。すごく痛いという訳ではないのだが、ロキソニンの薬効は5時間くらいしか持たないので、夜中に1度は目が覚めてしまうのが辛いところ。

抜歯後の食事は全てカロリーメイトのゼリーで済ませた。

抜歯翌日(2018/05/12)

今朝は特に痛みがないと思ったが、そういえば夜中に一度痛みで目が覚めてロキソニンを飲んだのだった。時間を確認していなかったが明け方だったのかもしれない。

昨晩寝るまでの間は若干の出血を感じていたが、もう完全に止まっている。

抜歯した場所に触れないように注意しながら歯磨きをする。こういうときはタフトブラシが便利だ。前回歯を抜いた時はネオステリングリーンという洗口剤を処方されたのだが、今回はなぜか処方されなかったのでコンクールFといういつも使っている洗口剤でそっとうがいをした。この洗口剤も抜歯後に使用される場合があるもので、アルコールが入っていないので刺激が少ないが口の中はすっきりする。長期間の使用で着色しやすいという副作用があるらしい。

朝はそれほどでもなかったがやはり薬効が切れると多少痛みが出てきて、耐えられないほどではないものの不快感があるので今日は11時と19時にロキソニンを服用した。少し頬に腫れを感じる。

痛み止めさえ効いていれば全く痛みがなく、昨日抜歯したことすら忘れる。普通に食事できそうな気もしたが、固いものを食べるとうっかり抜歯した側で噛んでしまいそうなのでうどんやグラタンなどの柔らかいものを注意しながら食べた。

抜歯翌日というのは割と大丈夫なのだ。前回もそうだった気がする。

抜歯2日後(2018/05/13)

また痛みで目が覚めて明け方にロキソニンを飲んだ。昨日時刻を確認し忘れたので時計を見たところ3:50だった。10分くらいで痛みが引いて寝た気がする。

今日は朝から若干痛みというか不快感がある。前回抜歯したときはあまり気にならなかったが、今回は外観から分かる程度に腫れてきた。

今日は9時と19時にロキソニンを飲んだ。ロキソニンの影響なのか腹の具合が悪い。飲めば腹が痛いし飲まねば歯が痛いし割と辛い。やはり2〜3日後がピークなのか。

食事は薬さえ効いていれば結構普通にできるようになっているが、念のため長めに茹でたうどんを連投している。

抜歯3日後(2018/05/14)

やはり夜中に目が覚める。今回は2時頃だった。寝ている間にロキソニンを摂取するソリューションが求められる。

朝の違和感がないので昨日よりは多少マシ、一昨日と同等といった状況。頬はまだ腫れているが、触った時に軽い痛みを感じるエリアは狭まっており、多少は回復しているものと思われる。昨日がピークだったと思いたい。

今日は10時と15時にロキソニンを飲んだ。薬が切れるとやはり痛い。朝痛んできたのでロキソニン飲もうと思ってたら会議が思ったより長引いて会議どころではなかった。今考えると別に会議途中だろうがサッと薬飲めばよかった。

抜歯4日後(2018/05/15)

ロキソニンが切れると目が覚めるのであれば、目が覚める時間にロキソニンを飲めば早朝に起床できるのではないか。そう思って昨夜12時頃にロキソニンを飲んだら6時に目が覚めた。抜歯は生活習慣の改善に使えるのではないか。しかし痛いので結局7時にはロキソニンを飲んだ。

朝若干雑に歯を磨いたら少し出血してしまいビビる。まだ注意しないと。

昼を過ぎても痛みが出なかったのでコイツは行けるかと思ったけど3時を過ぎたらまた痛みだした。なんやかんやで前回も1週間くらいは痛かった気がする。

抜糸5日後(2018/05/16)

結局今日も寝る前、朝、3時頃と3回ロキソニンを飲んだ。激痛ではないが、ロキソニンあるなら飲むわという程度には痛い。薬さえ飲んでしまえば無痛。

明日抜糸する予定なのだが大丈夫なのか。よく考えると前回抜糸したときも(まだちょっと痛いけど大丈夫なのか…)と思っていた気はする。前回抜糸したとき、抜糸自体が割と痛かった気がするのであまり行きたくないのだが、虫歯に端を発する長かった歯医者通いも明日で一旦区切りと思うと感慨深いものがあるような気もするし、ないような気もする。

抜糸6日後(2018/05/17)

今朝も起きると歯(のあった場所)が痛い。寝ている間に歯ぎしりをしているような気もする。

抜糸もめっちゃ嫌なのでめっちゃ嫌だなと思いながら大学病院に向かった。若干早く着いてしまったところ若干早く案内をしてくれた。

処置を待つ間、隣の席の人の話が聞こえてきたが、月1回佐渡から通ってるとかなんとかですごく大変そうな雰囲気を感じた。

肝心の抜糸は「前のオペが長引いている」とかいう理由で、抜歯したのとは別の先生がやることになった。前回のときもそうだったのだが、抜歯した人が抜糸まで面倒をみる的な決まりになっているのだろうか。

抜糸はやっぱりちょっと痛かったが、前回よりは痛くなかった気がする。でも処置後にうがいしたら思ったより出血してて吐いた水が赤く染まっており若干引いた。わざとそうしてあるのかもしれないが、歯医者の診療台についてるうがい用の流し、血の色がすごく目立つ配色だから出血してるとビビる。処置の時間は1分もかかってないと思う。

まだ痛みが続いているので、ロキソニンを10回分処方してもらった。

大学病院で行う処置はこれで完了したということで、もともと紹介状を出してくれた行きつけの歯医者に治療終わりました連絡をするのでちょっとお値段かかりますという話だった。多少余分にかかった気はする。連絡も業務なのでしっかり金とってって欲しい。

ともかくこれで歯医者通いは終わりなので、何事もなければ次回は定期検診だ。めでたい。

抜歯7日後(2018/05/18)

一週間経った。

糸を抜いて大分良くなったなという気持ちではあったものの、今日もロキソニンは7時、12時、20時と3回飲んだ。

ただ、歯を磨いたら血が出るとか痛いとかそういうことはもうない。昨日くらいまでは歯肉がフワフワしていたような気がするが、大分しっかりしてきた。

薬飲んでれば痛くないし薬が切れると痛いという代わり映えのしない日々が続いていてボチボチ飽きてきた感じはある。

抜歯8日後(2018/05/19)

朝起きた時に割と痛かった。今日は結局朝と夜の2回ロキソニンを飲んだ。3回が2回になったので大分良くなっている印象がある。

もう食べ物は完全に自由に食べられるが、右も左も親知らずを抜いた穴が開いている状態なので、何をどう注意しても穴に食べ物が入ってしまう。そして後に抜いた方に食べ物が入ると少し痛む。前回もこんなに長い間痛みがあっただろうかと思ったが、前回は抜歯10日後くらいに実家に帰ったら割と辛いカレーが夕食に出てきて割と痛かったような気がした。まだ8日目なのでこんなもんだろう。

ところで左側の親知らずを抜歯してから今日で丁度1ヶ月だった。左側の方はというと痛みはもう完全に無いが、若干突っ張っているような違和感はある。

抜歯9日後(2018/05/20)

朝起きてもほとんど痛くなかった。朝食を普通に食べても大丈夫だったのでコイツはいけるだろと思っていたが昼食を食べた当たりで若干痛みのような違和感が出てしまい結局ロキソニンを飲んだ。今日は1回で済んだので、ぼちぼち薬も卒業できそう。

ところで右の抜歯跡に比べると左の抜歯跡は大分食べ物が入り込む事がなくなっている様に思う。1ヶ月も経つし埋まってきているのだろうか。

抜歯10日後(2018/05/21)

起床後少しの間痛いような気がしたが、しばらくそっとしていたら治まった。どうも寝ている間の歯ぎしりとかが原因のような気がする。

もうそんな痛くないしいけるやろとか油断してたら昼食後に痛くなってしまった。食べ物が入ると痛い感じがする。耐えられないほどではなかったが、別に我慢しても良いことも無いので結局痛み止めを飲んだ。

抜歯11日後(2018/05/22)

基本的には痛くないのだが、ふとした瞬間から痛みだすことがある。気になって舌で触ってしまうのが良くない気がする。一旦痛くなるとしばらく痛いので結局今日も1回ロキソニンを飲んだ。薬を飲んでしまえばその後痛むことはない。

抜歯12日後(2018/05/23)

変化が無いわけではないが、ゆっくりと回復していくだけなのでいい加減書くことが無い。とりあえず朝起きたら痛いとかそういうことはもう無い。

しかしやはりふと痛みだす時があって、今日も1回だけロキソニンを飲んだ。

抜歯13日後(2018/05/24)

電動歯ブラシが抜歯した跡に当たるとさすがに少し痛い、というのが昨日までの状況だったが、今朝はそれもあまりなかった。

抜歯跡に食べ物が入ると多少痛いが、すぐにおさまるのでやっと痛み止めを飲まずに1日過ごせた。2週間も経過すれば平気になるようだ。

抜歯14日後(2018/05/25)

もはや食べ物が抜歯跡に入らなければ全く痛くないのだが、抜歯した方で普通に噛むという行為をしばらくしていなかったのでうっかり頬の内側を噛んでしまって少し痛かった。

抜歯15日後(2018/05/26)

もう食べ物が抜歯跡に入っても痛くない。そしてこうなってしまうともう書くことがない。今となっては食べ物が詰まりやすい穴が奥歯の奥に空いているだけだ。前述のようにこの穴が埋まるのには数ヶ月とかかかるらしいので、記録はこのへんで終わりにしようと思う。

2018/06/18

左の奥歯を抜いて約70日が経過したが、だいぶ穴が歯肉で埋まってきているようで、舌先で触ると分かる程度になってきた(最初はまたなにか食べ物が詰まっているのかと思った)。右側の抜歯跡が若干知覚過敏のようになっているらしく、冷たいものを飲み食いするとしみるのが気になる。

Rule of Three

調べ物をしていて、このような法則があることをふと思い出した。

病害虫の発生率がp%以下であるということを95%の信頼度で言うためには、3/p個の対象を調べて病害虫が存在しないことを示せば良い。

この法則はRule of Threeあるいは3の法則と呼ばれている*1

Rule of Threeは医薬品の稀な副作用を調査するという文脈の上で、次のように言われることもある。おそらくこちらの表現の方が有名だろう。

100人調べて1度も観測されなくても別の100人中3人に事象が生起する可能性があり、100 万人調べて0でもなお他の100万人中3人に生起する可能性は否定できない

導出

なんとなく「ポアソン分布から導出できる」みたいなアバウトな覚え方をしていて、いざやってみたら意外と手こずってしまったのでメモしておく。

まず前提条件を整理すると次のようになる。

  • 事象が起こる確率pはすごく小さい(と予想している)
  • 試行回数nは大きい
  • n回の調査で1回も事象は発生しなかった

pは二項分布に従うと考えられるが、pが大きくnが小さい二項分布はnp=\lambdaをパラメータとするポアソン分布で近似できる。ポアソン分布において事象がk回起こる確率は

P(k) = \frac{\lambda^k e^{-\lambda}}{k!}

だが、今回k=0なので、

P(0) = e^{-\lambda}

となる。ここで、\lambdaの上限を考えてみよう。\lambdaが大きくなるほどPは小さくなるだけで\lambdaはいくらでも大きくできることはできる。ただ、あまり小さなPを想定するのは不自然なので、P=0.05に対応する\lambdaあたりを上限としよう。

わざとらしい感じになってしまったが要するに\lambdaの片側95%信頼区間を求めるということだ。

0.05 = e^{-\lambda}

\lambdaについて計算すると、

\lambda = -\ln{0.05} \approx 3

となる。要するに、n回の試行で事象の発生回数が0だったとしても、\lambda=n回の試行における事象の発生回数の期待値の95%信頼区間はほぼ0〜3ということになる。

これが

100 万人調べて0でもなお他の100万人中3人に生起する可能性は否定できない

に対応する。ここでさらにnp=\lambdaとしてnについて整理すると

n = \frac{3}{p}

となり、これが

病害虫の発生率がp%以下であるということを95%の信頼度で言うためには、3/p個の対象を調べて病害虫が存在しないことを示せば良い。

に対応する。

参考

*1:私はこの法則が植物の病害虫調査で実際に利用されている例は見たことがない。住んでいた分野が違うだけかもしれないが

連続した空行(空白文字を含む)を置換する正規表現

連続する空行の置換

O'Reilly Japan - 詳説 正規表現 第3版に、連続した空行(これは、任意の空白文字を含む場合がある)を<p>に置換する正規表現として、

$text =~ s/^\s*$/<p>/mg;

というものが掲載されている(p.67)。これは、空白文字を含む連続した空行を含む文字列...with.\n\n\n\t \nTherefore......with.\n<p>\nTherefore...のように置換すると説明されていて、実際この例ならうまくいく。

3つ以上連続する改行が正しく置換できない

しかし、a\n\n\nbのように改行が3つ以上続く例はa\n<p><p>\nbと置換される。改行がいくつでも<p>は2つだ。『詳説 正規表現』は10年近く前の本なので、Perlの仕様変更があったのかもしれないけど、ちょっと調べた程度ではどうしてこのような挙動になるのかはよくわからなかった。ちなみにa\n\n\n\t\n\nba\n<p><p>\nbに置換されたりしてますますよく分からない。

パターン修飾子\mを付けた場合に連続する改行がどう解釈されるのか曖昧とかそんなことが理由じゃないかと思っていたけど、どうやらこれは/gを付けてグローバルマッチをさせた場合に起こるらしかった。/gを付けた場合に実際に何が起こるか、という点を正確に理解できていないのかもしれない。

/gを指定せずs/^\s*$/<p>/mとすれば予想通りの動きをする(もちろん1箇所しか置換できない)。

ということは、

while($text =~ s/^\s*$/<p>/m){}

とかやってやれば期待通りの置換ができる。

別解

とはいえ/gがあったりなかったりで挙動が変わるものを使うのは気持ちが悪いので、他の方法を探していたところStack Overflowに答えを見つけた(cf. Perl Regex To Condense Multiple Line Breaks - Stack Overflow)。

Perl 5.10から導入された文字クラスとして、\R\hがある。\Rは改行っぽいものにマッチする文字クラス(cf.perlrebackslash - Perl 正規表現逆スラッシュシーケンスとエスケープ - perldoc.jp)で、\hは水平方向の空白文字(タブとかスペースとか)にマッチする文字クラス(cf. perlrecharclass - Perl 正規表現文字クラス - perldoc.jp)であって、これを組み合わせると、連続した空行の置換は、

$text =~ s/\R(\h*\R)+/\n<p>\n/mg;

と書ける。\R\hは使えない言語もあるかもしれないが、こちらの方が何をやっているのか明確という印象はある。

他の言語

ところで、他の言語でも同じような現象が起こるかというとそうではないっぽい。途中で面倒になったのであんまり試してないけど。。

R

例えばR、gsub()で試してみると問題なく置換できる。

Python3

PythonもOK。

Go

Goもいけた。

JavaScript

試した範囲でPerlと同じ現象が起こるのはJavaScriptだった。/gをオプションで指定するタイプだと発生するっぽい?