データ分析メモと北欧生活

旧Untitled Note. データ分析、計量経済・統計とR、水産管理、英語勉強、海外生活などについて備忘録や自分の勉強のOutputの場所として

MENU

【ノルウェー】スタヴァンゲル旅行記2:スタヴァンゲル市街

スタヴァンゲル旅行記の続きです。


1つ目のプレーケストーレンについては以下で書きました。
keita43a.hatenablog.com


スタヴァンゲル市街についての動画も作ったので、ぜひご覧ください。

Stavanger スタヴァンゲル


スタヴァンゲルは1960年代の北海の石油開発で大きく発展した街で、その前は北海のニシン漁を中心とした漁業の町だったので、その頃の古い町並みと発展した新しい街が入り混じっているところがあります。

Gamle Stavanger (スタヴァンゲル旧市街)

1700~1800年代に造られた家が未だ残る歴史地区です。

木造で白い壁のかわいい家が立ち並ぶエリアなのですが、カラフルな花を軒先に植えている家が多く、とても映えていました。

f:id:keita43a:20200719064319j:plain

歴史地区ですが、出入りは自由です。
エリアの中にはギャラリーや、ガラス細工などのクラフトショップがあります。
今でも住んでいる方々がいるので、写真を撮るときはご注意を。

ちなみにこのエリアの中にあるノルウェー缶詰博物館に行きたかったのですが、コロナ禍の間にリノベーションしているらしく、2021年まで閉館していました。英語でのレビューでは、「石油博物館もよいけど、こちらのほうが面白い!」というコメントもあったほどなので、残念。

Skagen スカーゲン・ウォーターフロント

湾を挟んで旧市街と逆側のエリアは対象的に、古い建物を残しつつも近代的なウォーターフロントエリアになっていて、
レストランやバーが並んでいます。

f:id:keita43a:20200719065132p:plain

夏場はテラス席も多く出ているようですが、お店によってはヒーターを用意しているところもあったので、ノルウェーの夏の寒さに慣れない日本人でも安心です。

Fargegaten (Øvre Holmegate) ファルゲガーテン

Skagenより少し北東向きに進むと、地図上ではØvre Holmegateという名前のストリートがあります。
このエリアはFargegaten、ノルウェー語でカラフルストリートという愛称で呼ばれていて、
文字通りビビッドな色に染められた家々が立ち並ぶカフェ・バーエリアになっています。

f:id:keita43a:20200719065735j:plain

ノルウェー石油博物館

f:id:keita43a:20200719070929p:plain

ノルウェーは1960年代からの北海油田開発で大きく成長しましたが、スタヴァンゲルはその前線基地として発展してきました。
そんなノルウェーの石油開発の歴史が学べる博物館です。

単に石油の開発に関する工学的な解説だけではなく、石油で得た富をどう扱っているか?にも重きをおいています。
ノルウェーでは石油で得た資金を政府のファンドが運用していますが、運用先はノルウェー以外の海外のみとなっています。
これは、多くの資源国が陥った「資源の呪い」と言われる資源国の経済成長が止まる現象や「オランダ病」と言われる資源輸出に起因する為替レート上昇が製造業などの経済に影響する現象を避け、現世代が受けている恩恵を将来世代にも残すための知恵だとされています。

また、地球温暖化問題や、海底で作業するダイバーの健康問題など、石油開発の負の面にもしっかりと焦点をあてて解説しています。

f:id:keita43a:20200719070944p:plainf:id:keita43a:20200719071004p:plainf:id:keita43a:20200719071020p:plain

Byparken (ビーパルケン)

英語がわかる人だと「バイパーケン」と読みたくなりますが、ノルウェー語だとByはCityという意味なので、シティーパークという意味になります。

f:id:keita43a:20200719071219p:plain

中心にある湖の横にスタヴァンゲル大聖堂がありますが、今は工事中でした。

Sverd i Fjell (岩上の剣)

直訳すると山の上の剣、ですが、山というほどのところではないのでご安心を。

f:id:keita43a:20200719071705j:plain



スタヴァンゲル市街から少し離れたところ、車だと10分ぐらいで、バスだと20分ぐらいのところにあります。
歴史がありそうな感じですが、建てられたのは1983年。しかし、そのストーリーは北欧神話サガに基づいていて、ハラルド・フォールファグレによるノルウェーの統一が行われた872年の戦いがこのフィヨルドであったことを記念している銅像です。
3つの剣はそれぞれ、平和・融和・自由を表しています。

ノルウェーストーンヘンジ

f:id:keita43a:20200719211712j:plain

これはスタヴァンゲル市街ではなく、ヨルペラン(Jørpeland)というスタヴァンゲルとプレーケストーレンの間にある街にあるのですが、小さな島全体がJørpelandsholmenという公園になっているところがあり、たくさんのヤギが歩いています。

f:id:keita43a:20200719212805p:plain

もう一つの島の上に立てられた不思議なモニュメントがあるのですが、その形から「ノルウェーストーンヘンジ (Norwegian Stonehenge)」と呼ばれています。正式な名前はSolspeiletで、グーグルマップ上では英語でSun Mirror at Klungholmenとなっています。

f:id:keita43a:20200719212716j:plain

行き方が少しややこしくてGoogle Mapsでナビするとこの付近でバグってしまうので、スーパー(Rema 1000)を目的地にするとよいです。近所のスーパーやショッピングモールの駐車場に車を停めてから、徒歩で以下のボートセンターみたいなところの横の路地を抜けて橋に向かいます。

その他

他にも植物園にも行きましたが、ちょっと広めの公園という感じで、観光に来てわざわざ行くほどのところではなかったです。

今回は時間や天候の都合で行けていませんが、少し足を伸ばしたところにあるソラ・ビーチや、石器時代のファームバイキングハウスなんてのもあります。

すこしお金に余裕がある方には、庭園とディナー付きフィヨルドツアーなんていかがでしょうか?

【ノルウェー】スタヴァンゲル旅行記1:プレイケストーレン

ノルウェーではCOVID-19の感染もほぼ終息状態にあり、ソーシャルディスタンスや消毒などに注意しつつもほぼ普通の生活が戻ってきている。一応7月15日からヨーロッパ諸国に対して国境が開かれるのだが、さすがにこの時期に国外に行くのは色々リスクもありそうなので控えて、夏はノルウェー国内を旅行することにした。

というわけで、ノルウェー第4の都市スタヴァンゲル(Stavanger)に行ってきました。*1

ベルゲンから車でスタヴァンゲルに行って4泊という予定で行った。

行ったのは

  • プレーケストーレン (Preikestolen)
  • モナフォッセン(モナ滝, Månafossen)
  • スタヴァンゲル市街(オイルミュージアムなど)
  • フルォーリ(Flørli)

本当は崖の間に岩が挟まっているので有名なシェラーグにも行きたかったのだが、天候が悪く断念。しかし、全体として楽しい旅行だったので旅行記をつけておく。

プレーケストーレン

プレーケストーレンとは、リーセフィヨルドに面した大きな崖で、まさに断崖絶壁。
ノルウェーフィヨルド、というイメージ画像でもたびたび出てくる「ザ・フィヨルド」といえる場所である。


User:Smtunli - Wikimedia Commons


この上からの写真を撮るのはかなり危ないところを通るので、断念した。

アクセス

スタヴァンゲルからフィヨルド観光を含んだツアーや、バスが出ている。

私は車で行ったが、駐車料金が250NOK(約2800円)もした。

トレイル情報

Visit Norwayによれば所要時間は4時間とあるが、トレイル慣れしていない人はもう少し時間を見たほうがいいと思う。私たちはゆっくり上り下りして4時間ちょっとぐらいだった。英語サイトのOutttによれば2〜4時間とあるが、2時間は明らかにノルウェー人の時間なので当てにならない。*2

それなりに高低差があること(440m)に加えて、岩の多いトレイルで足場が悪いため、しっかりした登山靴を用意したほうが安全だ。

トレイルは、急峻な部分が2箇所ある。駐車場から森に入ってしばらく行くと、岩でできた階段のようなトレイルが続く。
f:id:keita43a:20200717042215p:plain
f:id:keita43a:20200717042233p:plain

しばらくすると湿地帯のようなところに出る。

f:id:keita43a:20200717042246p:plain

この湿地帯のあとにまた岩の階段のようなトレイルを登っていくと、大きな岩の上にようなところに出る。
ここからはほぼ平坦だが、岩の上を歩くような状態なので、天候によっては滑らないように注意が必要だ。
f:id:keita43a:20200717042306p:plain
f:id:keita43a:20200717042319p:plain
f:id:keita43a:20200717042425p:plain

しばらく歩いていると、フィヨルドに突き当たって到着。

f:id:keita43a:20200717042406p:plain

本当に崖の上で、フェンスも何もないので乗り出さない方がよい。
この写真はスマホだけ突き出して撮った。
f:id:keita43a:20200717042335p:plain

美しいリーセフィヨルドの景色が見られる。
f:id:keita43a:20200717042353p:plain


この写真は後日にフィヨルド観光船から撮ったものだが、まるで人工物のように角ばって飛び出している。


トレイルとしては人気で、子供連れや犬を連れた人も多く見かけた。
迷うような道でもないので、きちんとした靴と服装を準備すれば難しいハイキングではないのでおすすめである。

プレーケストーレンの動画も作ったので、よければご覧ください。

Preikestolen, Norway プレーケストーレン・ノルウェー

*1:都市圏としてはベルゲンに次ぐ第3位とされることもあるが、都市単体の人口はトロンハイムのほうが多い。

*2:ノルウェー人時間という概念はないが、ノルウェー人は一般的に早い。平均的に鍛えられている上に、スポーツ的にハイキングする人たちが一定数いるため、景色など楽しまず早さを意識する人もいる。

【ノルウェー】SbankenでBankIDが発行されない問題(銀行口座開設)

昔書いたこの記事の続きを書き忘れていたことに気がついた。

keita43a.hatenablog.com

問題:本人認証できないので口座開けない。

Sbankenで銀行口座を開こうとすると、BankIDが取得できない。(その後、口座すら開けなくなった)

ノルウェーでの公的なサービスを利用するためにはセキュリティのためにいくつかのIDが存在する。これは個人識別ID(Fødselsnummer)と紐付いているが、別物だ。

自分の税金の状態にアクセスする、ノルウェーの病院に予約を取ったりメッセージを送るHelsenorgeのアカウントにアクセスする、などのためにはセキュリティレベルの高いIDが必要で、その一つがBankIDである。BankIDは銀行口座を開けると取得でき、BankIDは一人に一つなので、いくつも口座を開いても複数のBankIDを取得できるわけではない。逆に、BankIDを持っていると別の銀行で口座を開くときに手続きが早い。

その他にもいろんなシーンで使うことがあるので、ノルウェーに住むなら必須のIDである。IDと言って物理的にカードなどが発行されるわけではなく、個人識別IDに紐付いたパスワードと銀行のワンタイムパスワード発行機で発行するパスワードを使うが、BankID取得後はスマホに登録して使えるBankID på Mobileがあるので、ワンタイムパスワード機を携帯する必要もなくなる。

店舗のないSbankenの問題。

流行りのネットバンクで手続きが早いと聞いていたのでSbankenで開いたのだが、問題は店舗のないことだった。実店舗がないので、本人の認証を郵便局(Posten)で行うというシステムになっている。アカウントを開くときに指示されたとおり、最寄りのPostenでIDをスキャンするのだが、当然IDとしては日本のパスポートしか無い。しかし、郵便局員がパスポートを何度やっても認識されないと言われる。

日本のパスポートはRFID入っているのに、と思ったけど、どうやらスキャン機自体が、ノルウェーで2010年4月以降に発行されたパスポートでないと認識されないらしい。

私が開いたとき(2018年秋)には、スキャン機が入れ替わった直後だったらしくて、郵便局員いわく”古い方法”で、物理的にパスポートのコピーを取って送付することで銀行口座自体は開けたが、それでもBankIDを取得するためにはスキャンが必要らしい。
また私の妻が後に開こうとした2020年の4月ごろは、その方法でも口座自体開くことができなくなっていた。

解決法:結局、実店舗のある銀行で開いた

ではどうするか?と、聞いたところ、公証人役場(ノルウェーでは裁判所?)に行き、パスポートが本人IDをであることを証明してもらう書類を作ってもらい、それを送れ、と言われた。

以下の記事でも、同様に苦労して公的書類を作った学生のエピソードが書かれている。
Exchange student struggles to open bank account - The Oslo Desk

すごくめんどくさそうだったので、結局DNBで作ってしまうことにした。
DNBならオンラインで申請したあと、実店舗でパスポートで認証できる。受け付けてからの時間はかかるが面倒が少ない。(BankIDを持っていれば手続きが早いらしい。BankIDを取得するための手続きなので仕方がない。。。)
留学生などはSpare Bank 1で作ってる人も多い印象。

他にどの銀行がいいのか?というのはこちらの記事も参考にご覧ください。

keita43a.hatenablog.com

【mac】Homebrewをインストールしたのに「command not found」が出る。

Homebrewをインストールしたのに、brewと打ってもコマンドが見つかりませんになったときの対処のメモ。

検索すると、「Homebrewをインストールしてないから」という回答が日本語でも英語でも出てくるのだが、そうじゃねえんだよと。

lanchesters.site


問題は、Homebrewコマンドラインでインストールできたはずなのに(brew helpと打ってみよう!と言っているのに)、いざ打ってみるとcommand not foundと出てくるということである。(スクリーンショット

f:id:keita43a:20200617211602p:plain


さて、注目すべきはWarning: /usr/local/bin is not in your PATH.という部分であった。
つまりパスが通ってなかったので、それを通してやればよかった。

keita$ export PATH=/usr/local/bin:$PATH

これでbrewが通るようになった。

参考:
github.com

【R】geom_sfで描画エラーが出たけどsfアプデしてないだけだった

簡単な備忘録。tidyverseやR本体などアップデートしたら、付随する他のパッケージもアップデートしましょう。

以前は普通に走ってたポリゴンファイルを読みこんでggplotで描画するコードにエラーが出たので調べたら、
tibble (というかtidyverse)はアップデートしたのにsfはアップデートされてなかったのでしましょう、というオチだった。


github.com

エラー

Error: All columns in a tibble must be vectors.
x Column `geometry` is a `sfc_GEOMETRY/sfc` object.
Run `rlang::last_error()` to see where the error occurred.

コード

jp_sh = readRDS("gadm36_JPN_1_sf.rds")

jp_sh2 = st_simplify(jp_sh, preserveTopology = TRUE, dTolerance = 0.01) 

ggplot(jp_sh2) + 
  geom_sf()

sfパッケージをアプデしたら普通にできました。

f:id:keita43a:20200607185220p:plain

他の可能性

他には、開発バージョンのtibble使ってたら同じエラー出たケースとかあるので、tibbleとsfの整合性が問題になる部分があるみたいです。

github.com

【R】新しい回帰分析表のパッケージ {modelsummary}

f:id:keita43a:20201209204033j:plain

Rで計量経済や統計分析やる時に、結果をきちんと理解しながらモデル作りたいし、できたモデルの結果を書き出すのも間違いなくやりたいですよね。

新しい回帰分析表のパッケージを発見したので、ざっと試してみました。

Rの回帰分析のパッケージ

Rでの回帰分析の表を書き出すのには {stargazer}というパッケージがあって、かなり柔軟にかつ自動的にhtmlやLatexに書き出しができるのですが、作者の都合か、更新が止まってしまってます。更新が止まってしまうと、新しい分析パッケージの結果とかを受け付けられなくなるんですよね。特に計量経済における回帰分析で役立つ{estimatr}の結果を受け付けないのがつらいです。

keita43a.hatenablog.com


他のパッケージで書き出すことも可能です。例えば{huxtable}パッケージのhuxreg関数は、かんたんに書き出せるし、Rのコンソールに見やすい結果を出力してくれるので良いです。しかし、まだ発展途上ということで、{stargazer}ほどの柔軟性はないのと、Latexで出力する際に表の位置が真ん中にいかないという問題がありまして、まだ実用的ではないかなという個人的な印象です。

以前なんとかestimatrの出力結果を固定効果付きでhuxregで書き出すという記事も書きました。

keita43a.hatenablog.com

今のところは作業中はhuxregで結果見ながら作業して、最終的に出力するときはいろんなワークアラウンドを駆使してstargazerできれいな表を書き出すという形をとっています。

{modelsummary}パッケージ

そんな時に、以下のTweetを見つけました。

新しい統計モデルのサマリー表ですって??
しかも、スレッドの次のTweetを見てみると、なかなかきれいな表が出力されています。


パッケージの説明はこちら。

vincentarelbundock.github.io


さっそくインストール

{modelsummary}がCRANに上がっているというツイートだったので、CRANからダウンロードしましょう。
dplyrの最新版に依存してるみたいなので、tidyverseも一緒に更新しておきましょう。

install.package("modelsummary")
install.package("tidyverse")

シンプルに回帰分析の出力

とりあえずmtcarsデータ使って、回帰分析を行って、結果を出力してみます。
結果はmsummary()関数を使います。

reg <- lm(mpg ~ cyl+disp, mtcars)

msummary(reg)

この関数を使うと、コンソールではなくViewerにきれいなテーブルが出力されます。

f:id:keita43a:20200529203525p:plain


2つ以上のモデルを同時に表示させる場合には、summary()のようにオブジェクトをそのまま並べるのではなく、リストにします。
リストの要素指定時にモデル名をつけることができて、これが列名になるようです。

regs <- list()
regs[['Model A']] <- lm(mpg ~ cyl+disp, mtcars)
regs[['Model B']] <- lm(mpg ~ cyl+disp+hp, mtcars)

msummary(regs)

f:id:keita43a:20200529203629p:plain

出力フォーマット

コンソールに表示

オプションでmarkdownと指定すれば、huxregのようにコンソールに表示することも可能です。

msummary(regs, 'markdown')
|            |Model A |Model B |
|:-----------|:-------|:-------|
|(Intercept) |34.661  |34.185  |
|            |(2.547) |(2.591) |
|cyl         |-1.587  |-1.227  |
|            |(0.712) |(0.797) |
|disp        |-0.021  |-0.019  |
|            |(0.010) |(0.010) |
|hp          |        |-0.015  |
|            |        |(0.015) |
|Num.Obs.    |32      |32      |
|R2          |0.760   |0.768   |
|Adj.R2      |0.743   |0.743   |
|AIC         |167.1   |168.0   |
|BIC         |173.0   |175.3   |
|Log.Lik.    |-79.573 |-79.009 |

出力フォーマットは以下のものに対応しているようです。

フォーマット オプション
HTML 'html'
LaTeX 'latex'
Text / Markdown / ASCII 'markdown'
RTF (Microsoft Word互換) 'ファイル名.rtf'
Images: JPG and PNG 'ファイル名.png'
gt テーブルオブジェクト(R) 'gt

指定もかんたんで、例えばLatexコードを出力したいときは

msummary(regs, 'latex')

と打つだけでOK.
ファイルに保存したい場合は

msummary(regs, 'result.tex')  # Tex ファイル
msummary(regs, 'result.png') # PNGファイル 

情報の編集

テーブル出力パッケージとして簡単にきれいな表を出力できることはわかったが、問題はどの程度柔軟に情報が編集できるかである。

見た目

まずテーブルの見た目については、modelsummaryはkableExtraやgtの柔軟性を応用して編集できるようになっているようだ。
私はgtはまだ使ったことがないが、テーブルを画像出力する時に力を発揮するらしい。
Rmarkdown/knitrで出力するケースではknitr::kableのエクステンションであるkableExtraが良い。
HTML用にもLatex用にも編集がしやすく、きれいなテーブルを出力できる。

htmlについてはgtパッケージが、latexについてはhtmlパッケージがデフォルトになっているらしいが、次のようなオプションで変更することも可能らしい。



例えばhtmlのテーブルをgtパッケージの関数で修飾してみると以下のようになる。

msummary(regs, output = 'html')%>%
  gt::tab_style(style = cell_fill(color = "lightcyan"), # 行に色を付ける
            locations = cells_body(rows = 3)) %>%
  gt::tab_style(style = cell_text(weight = "bold"),   #文字を太くする
            locations = cells_body(rows = 3))

f:id:keita43a:20200529203838p:plain

回帰分析表の要素

回帰分析表の要素については基本的な編集機能が実装されているようだ。

標準誤差・t統計量・p値・信頼区間

回帰分析表の係数の推定値の下に表示される統計量。経済学では標準誤差を示すのが一般的だが、他の統計量もstatisticオプションを変更することで示すことができる。

msummary(regs, statistic = 'std.error') # 標準誤差 (デフォルト)
msummary(regs, statistic = 'p.value')   # p値
msummary(regs, statistic = 'statistic') # t統計量
msummary(regs, statistic = 'conf.int', conf_level = .99) # 99% 信頼区間

f:id:keita43a:20200529203909p:plain
99%信頼区間を表示した場合

表のタイトルとメモ

表のタイトルもかんたん。titleオプションで指定します。

メモも簡単にnotesオプションで追加できるし、複数のメモをリストで表示できる。

msummary(regs,
         title = "表1:mtcarsの結果",          # タイトル
         notes = list('メモ1. 南無妙法蓮華経', # メモ
                      'メモ2. アーメン'))

f:id:keita43a:20200529203948p:plain

変数名の変更

変数名の変更は、coef_mapで行う。
予め変更前後の変数名の対応表をベクトルにしておくと便利そうだ。
新しく入るほうが右に来るらしい。(dplyrのrenameと逆なので注意)

# 変数名対応のベクトル
var_nam = c("mpg" = "Mile per gallon", "(Intercept)" = "Constant","cyl" = "Cylinder", "disp" = "Displacement",
            "hp" = "Horse Power")

msummary(regs, coef_map = var_nam)

f:id:keita43a:20200529195311p:plain

モデルフィットなどの統計量

msummary()ではデフォルトでR二乗などのモデルフィットの統計量が表示されている。
不必要な統計量はgof_omitオプションで非表示にできる。

msummary(regs, gof_omit = "R2|AIC|BIC")

f:id:keita43a:20200529200117p:plain

統計量のリストはgof_mapで見ることができる。

> modelsummary::gof_map
# A tibble: 14 x 4
   raw           clean         fmt   omit 
   <chr>         <chr>         <chr> <lgl>
 1 nobs          Num.Obs.      %.0f  FALSE
 2 r.squared     R2            %.3f  FALSE
 3 adj.r.squared Adj.R2        %.3f  FALSE
 4 AIC           AIC           %.1f  FALSE
 5 BIC           BIC           %.1f  FALSE
 6 logLik        Log.Lik.      %.3f  FALSE
 7 deviance      Deviance      %.2f  TRUE 
 8 df.residual   DF Resid      %.0f  TRUE 
 9 df.null       DF Null       %.0f  TRUE 
10 sigma         Sigma         %.3f  TRUE 
11 statistic     Statistics    %.3f  TRUE 
12 p.value       p             %.3f  TRUE 
13 df            DF            %.0f  TRUE 
14 null.deviance Deviance Null %.2f  TRUE
P値の星

stargazerと違って、msummaryではp値の星はデフォルトではつかない。(最近は星つけないジャーナルも多い)
starオプションをTRUEにすると、デフォルト(10%, 5%, 1%)の基準で星をつけてくれる。
もちろんカスタムも可能。"*"を変更することで好きな文字にすることもできる。

msummary(regs, stars = TRUE) # デフォルト, `* p < 0.1, ** p < 0.05, *** p < 0.01`
msummary(regs, stars = c("*" = .05, "**" = .01, "***" = .001)) # カスタム

f:id:keita43a:20200529200933p:plain
カスタムの結果

小数点以下の調整

小数点以下の調整はfmtオプションで行う。sprintfのフォーマットをとるので、
例えば小数点第2位までを表示したければ以下のように書く。

msummary(regs, fmt = '%.2f') # 小数点以下第2位まで


f:id:keita43a:20200529201440p:plain

huxregと違って、observationの数はきちんと整数になっているのがいいですね。

ちなみに、科学的表記法にする場合はfをeに変える。

> # 例
> sprintf("%.3f", 0.5342)
[1] "0.534"
> sprintf("%.3e", 0.5342)
[1] "5.342e-01
行の追加

例えば固定効果を入れたかどうか(Yes/No)などを書いた行を追加したい場合、stargazerではそれ専用のオプションが合ったが、
msummaryではhuxregのようにカスタムされた行を追加することで対応するらしい。
しかし、huxregより追加は容易である。

add_rowsオプションで追加する。表と同じ行数の文字ベクトルを作成して、listにして追加する。
(2020年9月19日追記) Ver. 0.4.1 より add_rows はデータフレームしか受け付けなくなったので、文字ベクトルではなくデータフレームで情報を与える必要があります。

##row1 <- c('固定効果1', 'No', 'No')  # Ver 0.4.1.以前
##row2 <- c('固定効果2', "No","No") # Ver 0.4.1.以前

# 0.4.1.よりデータフレームで与える
# tibble::tribbleで行ごとにデータを与えてtibbleを作れる。

rows <- tribble(~term, ~`Model A`, ~`Model B`,
               'FE 1', 'No', 'No',
               'FE 2', "No","No")  # Sep 20, added. Now the add_rows accepts only data.frame

attr(rows, 'position') <- c(15,16) # position attributeでどこに挿入するか指定する。15-16行目に追加

## msummary(regs, add_rows = list(row1, row2))# Ver 0.4.1.以前
msummary(regs, add_rows = rows)

f:id:keita43a:20200529202212p:plain

また、上の表な手作業ではなく、モデルをカスタマイズして固定効果を表示する方法も提案されていますが、別記事でまとめたいと思います。

で、estimatrは使えるのか?

大変使いやすそうだし、broom::tidyの出力に基づいているので柔軟そうだ。

ではestimatr::lm_robustは使えるのか?
結論から言うと、使えます!

regs_rob <- list()
regs_rob[["Model A rob"]] <- lm_robust(mpg ~ cyl+disp, 
                                       se_type = "stata",
                                       fixed_effects = ~ vs + am,mtcars)
regs_rob[["Model B rob"]] <- lm_robust(mpg ~ cyl+disp+hp, 
                                       se_type = "stata",
                                       fixed_effects = ~ vs + am,mtcars)

msummary(regs_rob)

f:id:keita43a:20200529202751p:plain


自動で標準誤差のタイプまで表示してくれていますね。

まとめ

  • {modelsummary}はシンプルで使いやすい回帰分析表出力パッケージ
  • broom::tidyの結果を受け付けるので、最近開発された統計モデルパッケージ等にも対応
  • gtやkableExtraを使ってテーブルが修飾できるし、様々なフォーマットへの出力に対応 (html, latex, RTF, text, 画像ファイル)
  • 多くのオプションでテーブルの内容もカスタマイズ可能

まだ実際に論文書くのには使ってはいませんが、よさそうなので使ってみます!
使ったコードはここにおいてます。

github.com

【R】スクレイピング時のループで404エラーが出た時に止まるのを避ける

ちょっと必要になって、ウェブサイト上のテーブルをまとめてダウンロードするという作業が必要になった。

だいたいPythonとかでやるのが主流なようだが、私はR使いなのでRでやる方法を探しているとrvestというパッケージがよいらしい。

rvestの基本的な使い方に関しては、Shohei Doiさんのノートが大変参考になった。
shohei-doi.github.io


さて、問題は同じようなウェブサイトページにURLをループで取得したい時に、時々ページ自体が存在せず、404 Error Not Foundのエラーが出ることがある。

これが困るのは、そのエラーが出た時点でループが止まってしまうことだ。とりあえず404があれば無視して次に行ってほしいので、どうすればいいかを調べてみたことをメモっておく。

tryCatchを使う。

以下のコードが実際に私が書いたコードの一部である。

ある記事に番号が振られており、1から300まで回して中にある表を取得しようとしているのだが、
不規則に404 Not Foundがでてくるので、単純にループを回すとそこで止まってしまう。

まず、delayedAssignでdo.nextという変数にnextの関数を割り振る。nextはloop内で出てくると、次のループに飛ばしてくれる命令である。

tryCatchは、1つ目の引数を評価して、エラーなど異常な条件が出ると、それに対して対応を決められる関数である。ここではrvestのread_htmlの関数を評価しているので、これにエラーがでると次のerror =で対応を決めている。
ここでは、もしURLが404 Not Foundでエラーが出た場合、"URL Not Found, skipping"というメッセージを表示して、force(do.next)により、次のループに飛ぶように命令している。

library(tidyverse)
library(rvest)

yr = 2017  #extracting year

list_tables = list()  # make list to store the result

for(i in 1:300){ # j-meldings up to 300 per year? maybe more in other years...
  
 url <- paste0("https://www.fiskeridir.no/Yrkesfiske/Regelverk-og-reguleringer/J-meldinger/Utgaatte-J-meldinger/J-",i,"-",yr)


 delayedAssign("do.next", {next})
 
 # use tryCatch to avoid stopping by error reading 404 page. 
 html <- tryCatch(
   read_html(url),
   error=function(e) {print("URL Not Found, skipping")
     force(do.next)})

 # check if <head> and <body> are properly read.
 
 # read html file from the url
 html <- url %>% read_html

 # print the number of meldings 
 print(paste0("melding ",i," is successfully read"))
 
# 2.1. check the title of melding
 
 # if the read html is not about the cod fishery, skip to next loop
 html_title <- html %>%
   html_node("h1") %>% 
   html_text() %>%
   str_detect(pattern = "Forskrift om regulering av fisket etter torsk, hyse og sei nord for 62°N")
 
 if(!html_title){
   next
 }