Yoshi Nishikawa Blog

医学となにかのインタラクティブ

Rを使って時系列データの変化点をみつける: changepoint

時系列データを扱う時に、変化点を見つけたいことがある。

Rにchangepointというパッケージがあるので、実装してみる。

  • 1つのchangepointはdefaultのmethod = AMOCで良い。

  • 複数のchangepointであれば、methodでPELT, SegNeighやBinSegを指定すれば良い。

#install.packages("changepoint")
library(changepoint)

# change in variance
x=c(rnorm(100,0,3),rnorm(100,0,10))
res1=cpt.var(x, )
plot(res1)
print(res1) 

# change in mean
y=c(rnorm(100,0,1),rnorm(100,7,1))
res2=cpt.mean(y)
plot(res2,cpt.col='blue')
print(res2)

# change in mean and variance
z=c(rnorm(100,0,1),rnorm(100,3,7))
res3=cpt.meanvar(z)
plot(res3,cpt.width=3)
print(res3)

# change in mean and variance, multiple changepoints
w=c(rnorm(100,0,1),rnorm(100,3,7), rnorm(100,-5,1))
res4=cpt.meanvar(w, Q=2, method = "PELT")
plot(res4)
print(res4)

f:id:yoshi_nishikawa:20171015153141p:plain

f:id:yoshi_nishikawa:20171015153200p:plain

f:id:yoshi_nishikawa:20171015153216p:plain

f:id:yoshi_nishikawa:20171015154151p:plain

計量経済分析ことはじめ~時系列データを正しく扱う~

Rによる計量経済分析を読んだ。このシリーズはとても勉強になる。

書籍のHPはこちら

キーワードと、キーとなるパッケージ・関数を列挙しておく。
時系列データを扱う時の問題点がわかりやすく解説されていて、よかった。

クロスセクションデータの回帰分析

最小二乗推定量 (OLS:ordinary least squares estimator)
決定係数 (coefficient of determination)
BLUE (the best linear unbiased estimator)

不均一分散

実世界のデータを扱う際に標準的仮定が満たされない場合がある。其の1つが不均一分散。

不均一分散の検定は以下で実行できる。

  1. Breusch-Pagan 検定

  2. White 検定

  3. Goldfeld-Quandt 検定: データを2分割して、それらを同じモデルを回帰した場合に誤差項の分散が等しくなるかどうかを検定

bptest(), gqtest()の説明書を見ながら、以下を実行。

library(lmtest)

## generate a regressor
x <- rep(c(-10,10), 100)

## generate heteroskedastic and homoskedastic disturbances
err1 <- rnorm(100, sd=rep(c(1,2), 50))  
err2 <- rnorm(100) 

## generate a linear relationship
y1 <- 1 + x + err1
y2 <- 1 + x + err2

## perform studentized Breusch-Pagan test
bptest(y1 ~ x)
bptest(y2 ~ x)

## perform Breusch-Pagan test
bptest(y1 ~ x, studentize = F) #studentize is a logical number
bptest(y2 ~ x, studentize = F)

## perform Goldfeld-Quandt test
gqtest(y1 ~ x, point=0.5) #point is interpreated to be a breakpoint.
gqtest(y2 ~ x, point=0.5)

時系列データの回帰分析

系列相関

時系列データでは、系列相関という問題が発生する。ある時点の回帰分析の誤差項と、其の隣の誤差項が相関している状態のことを指す。時系列データは、各時点毎に観測されたものであり、観測値の順番に意味があるので、注意が必要だ。

研究者は、このことをしっかり意識する必要がある。

系列相関の検定

  1. Durbin-Watson検定
  2. Breusch-Godfrey検定
library(lmtest)
## perform Durbin-Watson test
dwtest(y1 ~ x)

## perform Breusch-Godfrey test for first-order serial correlation:
bgtest(y1 ~ x)
## or for fourth-order serial correlation
bgtest(y1 ~ x, order = 4)

系列相関がある場合の推定

#FGLSを用いる方法
library(nlme)
res <- gls(y1~x, corr=corARMA(p=1))
summary(res)

#Newey-Westの修正を行う
library(sandwich)
coeftest(lm(y1~x), vcov=NeweyWest)

Stationary time series analysis (定常時系列分析)

定常過程 (Stationary process): 期待値・分散・共分散が時点tに依存せず、一定である。
independent and identically distributed; IID

モデル化には、自己回帰と移動平均を考慮して考える。

  1. Autoregressive model (自己回帰): p次のAR modelはAR(p)と記述される。

  2. Moving average model (移動平均

  3. Autoregressive moving average model (自己回帰移動平均

  4. Autoregressive integrated moving average model (自己回帰和分移動平均

ARCHとGARCH

時期によって分散の大きさがことなることを考慮する場合に用いる。

  1. (Generalized) autoregressive conditional heteroskedasticity model

GARCH(1,1)でモデリングしたいときは以下。

library(fGarch)
GARCH <- garchFit(formula = ~ garch(1, 1), data = data, 
         trace = FALSE, include.mean =TRUE)
summary(GARCH)

多変量時系列

複数の時系列データ(Aの株価とBの株価など)が互いに影響を及ぼしうるときに使う。

VAR vector autoregressive model

定常過程かどうか

library(tseries)
adf.test(data)
pp.test(data)

モデルの次数決定

library(vars)
VARselect(data, lag.max=5, type="const")
vars.res <- VAR(data, p=1, type="const")
summary(vars.res)
causality(vars.res, cause="A")
causality(vars.res, cause="B")

非定常時系列

単位根検定: ある時系列がランダム・ウォークに従っているかどうかの検定
見せかけの回帰: ランダム・ウォークに従っている者同士を回帰すると、見せかけの回帰(spurious regression)を来すことあり。 共和分分析: ランダム・ウォークに従っている者同士でも、共和分関係にあるものであれば、回帰にも意味がある。

library(tseries)
adf.test(data) #Augmented Dickey-Fuller
pp.test(data) #Phillips-Perron

library(urca)
ca.po(data) #Philips-Ouliaris
ca.jo(data) #Johansen

パネル

複数の同じ経済主体に関する複数時点でのデータである。

fixed effectとrandom effectの推定結果の検定: Hausman検定

library(plm)
phtest(res.fixed, res.random)

効率的なR運用を目指して

Rの基礎とプログラミング技法

Rの基礎とプログラミング技法に、Rを効率よく扱うコツについて非常にわかりやすく載っていたので、トライしてみた。

第5章:効率的なプログラミング

apply

  1. applyは、 ベクトル単位で処理できる
  2. lapplyはリストやデータフレームごと処理できる
  3. sapplyはlapplyの簡易版(返ってくるのはlistではなくベクトル)
  4. mapplyは多変量に適用できる
  5. tapplyはクロス集計表を作成できる
#apply

x<-matrix(c(1:9), nr=3)
apply(x,1,mean) #第二の引数は次元、ここでは行。2だと列。

#lapply, sapply

data(anscombe)
attach(anscombe)
ans.reg<-vector(4, mode="list")

for(i in 1:4){
x<-get(paste("x", i, sep=""))
y<-get(paste("y", i, sep=""))
ans.reg[[i]] <- lm(y~x)
}

lapply(ans.reg, coef) #データフレーム、リスト、ベクトルなどの構造に、coef()関数を適用
sapply(ans.reg, coef) #sはsimplifyであり、返されるのはベクトル

lapply(ans.reg, "[", 1) #このようにも[]を関数としても使える。

#mapply

mapply(sum, 1:10, 2:11, 3) #多変量に同じ処理が出来るが、一つ目の引数が関数となる。

#tapply

data(iris)
attach(iris)
tapply(Sepal.Width, Species, mean)

data(iris)
attach(iris)
tapply(Sepal.Length, Species, range)

効率性を分析する

以下のようなツールが搭載されている。

proc.time() #Rを起動してからの経過時間を出力する
system.time() #式の実行時間を出力する

Rprof() #コレ以降の実行時間を分析開始
Rprof(NULL) #これまでで終了
summaryRprof() #summary表示

gc() #R実行中のメモリの消費状況を調べるには、gc()を使う。其の点で使われていないメモリ領域が解放される。

さらに詳しい分析は、proftools packageなどが使える。

バッチ処理(第9章:拡張)

一括処理とも呼ばれる。定期的に実行しなければならない処理や、エンジニア同士でシェアするような状況では、Rを立ち上げずとも、以下のようにOSのバックグラウンドで出来るようにすべきである。

例えば、xxxx.rというファイルをRで実行して、output.txtを出力するには、以下をターミナル(macでは)で実行すればよい。

R CMD BATCH xxxx.r output.txt

分散分析についてANOVA(Analysis of Variance)

分散分析

数学いらずの医科統計学で、統計の概念的な理解を深めているのですが、 分散分析の理解については、以下の書籍や、リンクなどを参考に数式含めて再学習しました。

この東京大学出版会シリーズ、赤の統計学入門、緑の人文・社会科学の統計学と並んで、よく纏まっています。

あくまで、数式は使わずに言葉で纏めておきます。

一元配置分散分析

1つの要因への従属変数への影響を分析する(3グループ以上でも可能)。

http://elsur.jpn.org/resource/anova.pdf

分散分析1(一元配置分散分析)

私のための統計処理 ー基礎解説ーANOVA

  1. 平方和を分解
  2. 平均平方を算出
  3. 検討したい要因についてF(変動部分と誤差部分の二つの分散比)を求め、帰無仮説の棄却の有無を判断

これらのばらつきが有意に差があるなら、「平均値には差がある」という判断ができる。

例えば、A組、B組、C組の生徒に身長差があるかどうかを検証する。

多元配置分散分析

複数の 要因 (2 つ以上) による従属変数への影響を分析する。

この場合に、要因間の交互作用(interaction)を考慮する必要がある。

http://www.u.tsukuba.ac.jp/~hirai.akiyo.ft/forstudents/eigokyouikugaku7.files/2013_5_28.pdf

Clinical Pictureを投稿する

前回、レターを投稿するで好評いただいたので、clinical pictureに関しても纏めておきます。 *2019年10月19日Radiology追記しました。

投稿先あれこれ

NEJM

Please include a title for your submission. The title should contain no more than eight words. No more than two authors may be listed. Please provide the name, highest academic degree, address, e-mail address, telephone number, and fax number of each author. The legend should contain no more than 150 words.

http://www.nejm.org/page/author-center/images-in-clinical-medicine

JAMA

"What Would You Do Next?" with 4 single-phrase plausible treatment options describing possible courses of action with 1 being preferred Case presentation: 250 words Discussion: 500-600 words ≤10 references ≤3 authors 1-2 small figures Patient permission may be needed

http://jamanetwork.com/journals/jama/pages/instructions-for-authors/#SecArticleTypeClinicalChallenge01

Gastroenterology

Clinical Challenges and Images in GI (See full information about article types here) This article type is presented as an unknown with the diagnosis hinging on the correct interpretation and integration of the image and clinical data. Submissions must adhere to the following guidelines: Manuscript: only Microsoft Word documents will be accepted. Title Page: title (cannot reveal diagnosis); authors' names (limit of three); authors' institutions; corresponding author contact information; conflict of interest statement (for all authors). Word Count: Q&A format (one page each). Abstract: not required. Tables/Figures: no limit. Please submit figures as separate attachments in JPEG, TIFF, EPS, or PDF formats (300 PPI resolution). Each figure should be labeled with letters only (i.e., no subpanels). All necessary information pertaining to figures must be included in the text (i.e., no figure legends). The inclusion of three or more images will qualify the submission for online-only publication. References: limited to three.

Clinical Gastroenterology and Hepatology

Image of the Month Image of the Month presents a striking clinical image(s) that is meant to challenge and inform the reader. Although priority will be given to exceptionally unique submissions, and those that are not similar to recently published cases; authors should be encouraged to present quality images of more commonly encountered diseases and conditions. All submissions should contain no more than four (4) color or black and white images of 300 dpi resolution. Images should be accompanied by a Word document containing a brief description of no more than 200 words. The text should succinctly present relevant clinical information, including a short description of the patient's history, relevant findings, clinical course, response to treatment, and condition at last follow-up. No more than three authors are allowed on each submission. Contributors must provide their names, addresses, phone, and e-mail addresses.

http://www.cghjournal.org/content/authorinfo

Canadian Medical Association Journal

Clinical images Images are chosen because they are particularly intriguing, classic or dramatic. Preference is given to common presentations of important rare conditions and important unusual presentations of common problems. Structure: 300 word limit; up to two authors; up to three references formatted in the Vancouver style figure and case pertaining to a real patient; clear, appropriately labelled, high-resolution images must be accompanied by a figure caption A brief case description is followed by a concise explanation of the educational significance of the images that typically includes epidemiology, differential diagnosis, investigations, management and prognosis For examples, see Clinical images

Submission guidelines | CMAJ

Cleveland Clinic Journal of Medicine

http://www.mdedge.com/ccjm/page/clinical-picture

Netherland Journal of Medicine

Photo quiz A photo quiz should not exceed 500 words and include no more than two figures and four references conform the Vancouver style. Abbreviations of measurements should be quoted in SI units.

QJM(Quarterly Journal of Medicine)

Clinical pictures should be no more than 500 words. Ideally, we suggest 200-300 words. Only one figure file is allowed for the image. The selected figure file may contain two images (labelled figure 1a and 1b). Clinical pictures should have six references or fewer. The patient's consent is required for photographs of the patient. Consent should be in the form of a letter (written by author or patient) and must be signed by the patient.

Instructions to Authors | QJM: An International Journal of Medicine | Oxford Academic

Journal of General Internal Medicine

Clinical images report on visual findings in clinical medicine that have educational value. They can include radiology results, high quality clinical images, or electrocardiograms. Images should have a text description that does not exceed 200 words. No more than three authors may be listed. In the initial submission (for clinical images ONLY), each image should be sent as a separate file with the submitted text.

https://jgimed.org/authors/JGIM%20Instructions%20for%20Authors.pdf

European Journal of Internal Medicine

Internal Medicine Flashcards Please see Introducing the "Internal Medicine Flashcards": Call for papers Volume 24, Issue 6, published online only. Authors: maximum 3 Image: one, single or multi-paneled. Only original, high-quality images will be considered for publication, provided they do not contain material that has been submitted or published elsewhere. If a photo of an identifiable patient is used, a specific release form must be completed and signed by the patient and enclosed to the submission. All the printed information that might identify the patient or the authors' institution (including but not limited to the hospital or patient name, date or place) should be removed from the images Main section (case description): maximum 175 words Discussion section: maximum 225 words + maximum 3 references

Guide for authors - European Journal of Internal Medicine - ISSN 0953-6205

http://www.ejinme.com/article/S0953-6205(13)00206-9/fulltext

Internal Medicine

Pictures in Clinical Medicine Images depicting an instructive case, a novel finding and/or element of current interest for clinical medicine. The manuscript text should not exceed 150 words with up to 2 references. The picture legend should not be contained. The title should be concise and contain up to 8 words and there should be no more than 4 pictures. Only 4 authors may be listed.

Instructions to Authors|Internal Medicine

Radiology

仙台厚生病院 齋藤宏章 先生に教えていただきました。ありがとうございました。

Radiologyに掲載された、齋藤先生の論文はこちらです。

Images in Radiology
Images in Radiology are state-of-the-art images covering a range of both common and uncommon radiologic conditions. The imaging findings are compelling, visually appealing, and demonstrative. Images in Radiology captures the excitement of state-of-the-art imaging and discovery in our discipline. Images are of interest to a broad range of the readership. Images in Radiology are not intended to be case reports.
Original, high-quality images are considered for publication provided they do not contain material that has been submitted or published elsewhere. To submit an image for publication in Radiology, please follow the submission instructions below. At the discretion of the Editor, images may appear in the print version of the Journal, the electronic version, or both.
Submission
Title: no more than 8 words
Authors: no more than 2 authors
Figure Limit: 1-2 maximum
Legend word count: no more than 150 words
Figures and legend text should be in one double-spaced Word document
Figures must also be submitted in separate TIF, EPS, AI, or PSD files with 600 dpi
Include a Full Title Page

https://pubs.rsna.org/page/radiology/author-instructions

R markdownのknitr::opts_chunk設定をどうしているか。

忘備録的に。 私は、基本はこのスタイルで、

  1. figureの縦横比が一致するようにしている。
  2. tidyverseが静かに発動するようにしている。
knitr::opts_chunk$set(echo = TRUE, out.width = 480, out.height = 480, fig.width = 7,fig.height = 7)
library(tidyverse, quietly = T)

他にも、何か良いsetupあれば教えてください!