読者です 読者をやめる 読者になる 読者になる

りんごがでている

何か役に立つことを書きます

EzXML.jlを作った話

Julia

先月ですが、EzXML.jlというXMLやHTMLを扱うためのパッケージをリリースしました。 EzXML.jlはlibxml2というC言語で書かれたメジャーなライブラリのラッパーなのですが、それをJuliaでラップする過程で色々と工夫がありましたので、その辺を共有したいと思います。また、パッケージの開発についても私の思うところを色々と紹介していきたいとも居ますので、これからJuliaのパッケージを作ってみようと思う方には参考になるかと思います。

何故作ったか

そもそも、JuliaにはLightXML.jlというパッケージがあり、こちらもlibxml2をラップしたパッケージでEzXML.jlと重複する部分がかなりあります。LightXML.jlはJuliaのXMLを扱うパッケージの中では最もよく使われているもので、多くの依存パッケージがあります。しかしながら、いくつか問題があり、私はこのパッケージにパッチを送るでもなく新しいパッケージを作ることにしました。

LightXML.jlで問題だと思われたものを列挙すると主に以下の様なものが挙げられます。 * APIがJuliaっぽくない * 型の設計があまりうまくない * 機能が少ない * C側で確保されたメモリーを自動的に解放してくれない

まず、APIに関してですが妙な名前の関数が提供されています。例えば、XMLファイルをパースして読み込む関数はparse_fileという何をするのか分からない名前です。Juliaではしょっちゅうusingを使ってパッケージから提供されている関数をパッケージ名を付けずに利用しますから、例えばちょっと大きめのプログラムでparse_fileという関数が呼ばれていても、これがXMLを読み込む関数だとはなかなか思えません。また、幾つかの関数名にはchild_nodesis_elementnodeなどアンダースコアが使われていて、これもJuliaっぽい関数の命名法ではないです。例えば、Juliaの標準関数を見てみるとeachlinehaskeyなどアンダースコアを使わない関数名が一般的で、スタイルガイドにもそのように記述されています。

型に関しても、ちょっと設計がよくない気がします。XMLのDOMモデルでは、XMLの文書中を構成する要素をノードと呼び、<element/>のような要素ノードや<!-- comment -->のようなコメントノードなど様々な種類のノードがあるのですが、LightXML.jlではそのようなノードを別々のJuliaの型でラップしています。ですので、以下のようにノードを見つけるとその種類を判別して目的の型のノードのキャストする必要があります。

for c in child_nodes(xroot)
    if is_elementnode(c)
        e = XMLElement(c)
        # ...
    end
end

これはいかにも面倒くさいですし、一度生成したオブジェクトを別のオブジェクトへキャストすることで余計なオブジェクトの生成も起きますから、パフォーマンス上もあまりよろしくないと思います。

また、libxml2は名前空間XPathなど様々な機能を実装しているのですが、LightXML.jlから提供されているのはそのごく一部でXMLの読み書きとDOM操作だけです。それに、libxml2を通して確保したメモリー領域は、Julia側からの参照がなくなっても自動的には解放してくれないので、自分で解放しないとメモリーリークを起こします。

これらの問題を解決して快適にXMLを扱えるようにするため、EzXML.jlという新しいパッケージを開発することに決めました。

EzXML.jlの売り

EzXML.jlの売りは、その名前の通り簡単に(easy)にXMLを扱えるようにすることです。そのため、APIや型の設計をもう一度考え直し、かつlibxml2が提供する便利な機能を使えるようにしています。次のコードは、EzXML.jlを使ってXMLのツリーを辿ったり、XPathを使うデモです。

using EzXML

# Parse an XML string
# (use `readxml(<filename>)` to read a document from a file).
doc = parsexml("""
<primates>
    <genus name="Homo">
        <species name="sapiens">Human</species>
    </genus>
    <genus name="Pan">
        <species name="paniscus">Bonobo</species>
        <species name="troglodytes">Chimpanzee</species>
    </genus>
</primates>
""")

# Get the root element from `doc`.
primates = root(doc)

# Iterate over child elements.
for genus in eachelement(primates)
    # Get an attribute value by name.
    genus_name = genus["name"]
    println("- ", genus_name)
    for species in eachelement(genus)
        # Get the content within an element.
        species_name = content(species)
        println("  └ ", species["name"], " (", species_name, ")")
    end
end
println()

# Find texts using XPath query.
for species_name in content.(find(primates, "//species/text()"))
    println("- ", species_name)
end

EzXML.jlがXMLファイルを読み込むときの関数名はreadxmlです。これは、標準ライブラリにもあるreadcsvに沿って命名したものです。readという関数名はJuliaの標準ライブラリではファイルから読み込む際に使われる一般的な動詞ですので、その点でも一貫性があります。文字列からデータを読み込む場合にはparsexmlという関数名になっていて、これもJuliaではparseが文字列に対して頻繁に使われることを反映されています。ある要素の属性を取るには辞書(Dict)と同じようにelement[name]という構文が使えます。値をセットするにはもちろんelement[name] = valueと書きます。他の関数に関しても、一貫性があって一度理解すればあまりレファレンスを引かなくてもプログラミングできるようになっています。

型の設計に関してもなるべくシンプルになるようにしています。EzXML.jlが提供している型は主にEzXML.DocumentEzXML.Nodeの2種類だけです。EzXML.DocumentXML文書全体に対応する型で、EzXML.Nodeはその構成要素であるノードを表す型です。EzXML.NodeXML中のノードっぽいものすべてをラップしますので、要素も、テキストも、コメントも、属性もすべてがEzXML.Nodeです。そのノードがどのような種類のノードか知るにはnodetype関数を呼びます。これらの型の名前はコード中で使われることはあまりないですから、パッケージとしてはexportしていません。しかし、必要であれば上のようにモジュール名を付けて使うこともできます。

機能もLightXML.jlより豊富です。名前空間を扱うこともできますし、XPathを使って目的の要素を検索することもできます。XPathfindという関数をオーバーロードしていて、find(node, "//foo")とすればnodeノード以下のfoo要素をすべて取得することができます。また、メモリに乗り切らない大規模なXMLファイルを高速にパースするために、streaming readerも提供しています。

最後に、EzXML.jlはメモリリソースを自動的に解放してくれます。普通にJuliaを使っている人はメモリリソースの解放なんて気にしたくないでしょうから、自動的にやってくれれば大変楽です。なので、EzXML.jlではこの辺の面倒を見てあげることにしています。この次では、この実装がどのようになっているのかを説明していきます。

内部実装

まずはJuliaでCのライブラリーをラップする際の基本的なことを確認しておきましょう。JuliaのオブジェクトとC言語の構造体ではメモリレイアウトに互換性があります。例えば、struct point_s { int kind; double x; double y; }という構造体があれば、Julia側からはimmutable Point; kind::Cint; x::Float64; y::Float64; endのような型のオブジェクトを定義してやれば、そのメモリレイアウトをそのままJulia側に写すことができます。Cの関数からpoint_s*のようなポインタ値が返された場合は、Julia側でunsafe_load(ptr::Ptr{Point})を使ってJuliaのPointオブジェクトに変換できます。また、C言語の関数の呼び出しはccallで実現できます。

EzXML.jlでは、libxml2のノードを扱う時は以下の_Node型に写しています。これは、libxml2のどんな種類のノードも共通して持っているフィールドです。後で述べますが、最初の_privateフィールドが重要です。

immutable _Node
    _private::Ptr{Void}
    typ::Cint
    name::Cstring
    children::Ptr{_Node}
    last::Ptr{_Node}
    parent::Ptr{_Node}
    next::Ptr{_Node}
    prev::Ptr{_Node}
    doc::Ptr{_Node}
end

Juliaはオブジェクトが解放されるときに関数を呼び出す仕組み(finalizer関数)を提供しており、これを使うことでリソースの解放ができます。しかしこれは、想像するよりちょっと複雑です。何故なら、XMLのノードはlibxml2内で繋がっていますから、誤ってつながっている一部のノードを勝手に解放してしまうと、一貫性が保たれなくなってしまいます。以下のコードを見てみましょう。

function extract_first_child(filename)
    doc = readxml(filename)
    r = root(doc)
    c1 = firstchild(r)
    return c1
end

# foo.xml:
# <root>
#   <child1/>
#   <child2/>
# </root>
child = extract_first_child("foo.xml")

この関数ではXML文書中の最初の要素を返しています。関数が返った後はJuliaレベルではdocへの参照が無くなりますから、Juliaのガーベッジコレクション(GC)はdocを解放しようとします。このときにlibxml2のメモリーも解放して良さそうに思えますが、実際にはそうではありません。docが参照しているlibxml2の構造体はrを通してc1につながっていて、c1docへのポインターを持っています。ですので、領域を解放するには、XMLのツリーのどこにもJuliaからの参照が残っていない状態である必要があります。

f:id:bicycle1885:20161215084858p:plain

この問題を回避するために、EzXML.jlではメモリーの解放は常に特定のノードから起きるようにしています。実装上は、Node型のオブジェクトはptrというlibxml2のノードを表す構造体へのポインターと自身を管理しているownerという2つのフィールドを持っています。このownerオブジェクトが他のノードのメモリー管理を担っているオブジェクトです。

type Node
    ptr::Ptr{_Node}
    owner::Node
end

Nodeコンストラクターの内部では、以下の様なコードを使って最上位のNodeポインターを手繰り寄せて、そのノードのownerに指定しています。つまり、オーナーはXMLツリーの最上位にあるノードがその役割を果たします。

# determine the owner of this node
owner_ptr = ptr
while unsafe_load(owner_ptr).parent != C_NULL
    owner_ptr = unsafe_load(owner_ptr).parent
end

このownerフィールドを持つ理由は他にもあります。ownerを持つことで、GCがその子(子孫)ノードより先にオーナーを解放しようとするのを防ぐことができます。上の関数の例で言えば、c1docへの参照をJuliaのレベルで保持しているため、c1がある限りdocが先に解放されることはありません。さらに途中に挟まれているrは自身や他のノードのオーナーではないので、JuliaはrGCするときにlibxml2で確保されたメモリー領域を解放しません。こうして、常にノードは最上位のノードへの参照を持ち、最上位のノードがメモリーを解放することにして、誤ってGCされるのを防いでいます。

しかし実は、もうひとつ考えなければいけない問題があります。それは、Julia側からあるlibxml2のノードに対して複数の参照を作ってしまうことです。以下の様な簡単なコードを考えてみると、すぐに問題が分かります。

doc1 = readxml("foo.xml")
doc2 = document(root(doc1))

このとき、単純に実装するとdoc1doc2はlibxml2の同じノードを指しているのに別々のJuliaオブジェクトになっています。すると、あるノードの対して2つのオーナーが存在することになり、両者が最終的に同じノードを解放しようとして後のほうが不正な操作になります。

この問題はlibxml2のノードの構造体にJuliaオブジェクトへのポインターを保持することで解決しています。既に述べたように、libxml2のノードの構造体には_privateと呼ばれるフィールドがあり、ユーザーが自由に使うことができます。そこで、このフィールドにJuliaのNodeオブジェクトに対するポインターを保持しておき、あるノードに対するJulia側のNodeオブジェクトを作る際には、もしあればここからオブジェクトを取り出しています。つまり、既にJulia側であるノードに対するオブジェクトを作成していればそれを使い、無ければ新しく作ってポインターをノードの保持しておくという操作をNodeコンストラクター内で行っています。こうすることで、あるlibxml2のノードに対するJulia側のNodeオブジェクトは常に高々1個になり、二重に領域を解放してしまうことを防ぐことができます。

最終的に、Nodeコンストラクターは以下の様な感じになります(実際のコードを少し簡略化しています)。 CのポインターからJuliaのオブジェクトを取り出すのがunsafe_pointer_to_objrefで、Juliaのオブジェクトのポインターを取得するのがpointer_from_objrefです。

function Node(ptr::Ptr{_Node})
    # return a preallocated proxy object if any
    str = unsafe_load(ptr)
    if str._private != C_NULL
        # found a valid proxy
        return unsafe_pointer_to_objref(str._private)::Node
    end

    # find the owner of this node
    owner_ptr = ptr
    while unsafe_load(owner_ptr).parent != C_NULL
        owner_ptr = unsafe_load(owner_ptr).parent
    end

    if ptr == owner_ptr
        # manage itself
        node = new(ptr)
        node.owner = node
    else
        # delegate management to its owner
        owner = Node(owner_ptr)
        node = new(ptr, owner)
    end

    # register finalizer and store the pointer to a node object
    finalizer(node, finalize_node)
    unsafe_store!(convert(Ptr{UInt}, ptr), convert(UInt, pointer_from_objref(node)))

    return node
end

finalize_nodeは自身がオーナーである場合は、それ以下の管理するノードにつながっているJuliaのオブジェクトを切り離し、libxml2の関数を呼んでノードを解放しています。もしオーナーでなければ、_privateNULLポインターを代入して、JuliaのオブジェクトがGCされてもう存在しないことを示しています。

# Finalize a Node object.
function finalize_node(node)
    node_ptr = node.ptr
    if node === node.owner
        # detach pointers to C structs of descendant nodes
        traverse_tree(node_ptr) do ptr
            str = unsafe_load(ptr)
            if has_proxy(str)
                # detach!
                unsafe_extract_proxy(str).ptr = C_NULL
            end
        end
        # free the descendants
        if unsafe_load(node_ptr).typ == DOCUMENT_NODE
            ccall((:xmlFreeDoc, libxml2), Void, (Ptr{Void},), node_ptr)
        else
            ccall((:xmlFreeNode, libxml2), Void, (Ptr{Void},), node_ptr)
        end
    elseif node_ptr != C_NULL
        # indicate the proxy does not exit anymore
        store_proxy_pointer!(node, C_NULL)
    end
    return nothing
end

以上で、EzXML.jlの大まかな内部実装がお分かりになるかと思います。実際には、サブツリーを別のツリーへと動かす際にはオーナーをアップデートする必要があるなど細かい処理が多少あるのですが、そのへんまで興味ある方は元のコードを読んでみて下さい。まだ新しいパッケージですので、使ってみてバグや提案などがある場合は報告して下さい。

RユーザーのためのJulia100問100答

R Julia

R Adevnt Calendar 8日目の記事です。大幅に遅れて大変申し訳ないです。

この記事ではR言語ユーザーのために100問100答形式でJuliaを紹介していこうと思います。

Julia言語

Juliaってどういう言語なの?

Juliaは高レベルでハイパフォーマンスな技術計算のための動的言語だよ。書きやすさと実行速度の両立がウリの言語だよ。

誰が作ってるの?

主にボストンのMITの人達が中心に作っている言語だよ。特にJeff Bezonson, Stefan Karpinski, Viral Shah, Alan Edelmanの4人が初期の重要人物だよ。

自由に使えるの?

Juliaの処理系はMITライセンスで配布されているから、商用でもなんでもかなり自由に使えるよ。

どれくらい速いの?

すごく速いよ!大体C言語の2倍以内くらいの収まる速度だよ。

Rと比較してどうなの?

数倍から数百倍くらい高速かな。特にループや関数呼び出しがたくさんある場合には顕著な差が出るよ。

公式サイトにあるベンチマークの抜粋を載せておくよ (C言語=1.0)。 f:id:bicycle1885:20161214151939p:plain

実例はある?

forループがあるような場合はかなりはっきり差が出るよ。 例えば1からNまでのコラッツ予想を実証する関数collatzをJuliaとRで書いてみたよ。

function collatz(N)
    for n in 1:N
        while n != 1
            if iseven(n)
                n = div(n, 2)
            else
                n = 3n + 1
            end
        end
    end
end
collatz <- function (N) {
    for (n in 1:N) {
        while (n != 1) {
            if (n %% 2 == 0) {
                n = n %/% 2
            } else {
                n = 3 * n + 1
            }
        }
    }
}

N=100,000で試してみるとJuliaはRの300倍ぐらい高速だったよ。

julia> @time collatz(100000)
  0.026903 seconds (4 allocations: 160 bytes)

R> system.time(collatz(100000))
   user  system elapsed
  8.429   0.028   8.479

どうして速いの?

LLVMというコンパイラの基盤ライブラリをつかって実行時にコンパイルを行っているからだよ。あとJulia言語自体が最初からパフォーマンスを考慮して設計されているよ。詳しくはBezansonら論文を参照してね。

どこが技術計算に向いてるの?

実行速度が速く、かつ動的プログラミング言語で気軽に実行できるから短いスクリプトで計算するのに向いてるよ。あと、数値計算でよく用いられるデータ型や関数が最初から用意されているから、インストールしてすぐ使い始められるよ。例えば、線形代数のためにBLASLAPACKの関数が標準で入ってるよ。

Juliaのレポジトリはどこ?

GitHubのレポジトリがここにあるよ: https://github.com/JuliaLang/julia

インストール用のバイナリは?

公式サイトのダウンロードページにWindows, Mac, Linux用のが用意してあるよ: http://julialang.org/downloads/

Juliaのマニュアルやチュートリアルは?

Juliaはドキュメントがよく書かれているから、公式マニュアルを読むのがオススメだよ。 日本語だと私が書いたJuliaのチュートリアルM.Hiroiさんのサイトがあるよ。

環境

どうやって実行するの?

次のスクリプトをhello.jlとしてファイルに保存したとしよう。

println("hello, world")

Juliaをインストールしてあるならjuliaコマンドがあるはずなので、これにさっきのファイルを渡せばRscriptみたいに実行できるよ。

$ julia hello.jl
hello, world

対話セッション(REPL)を使うには、引数無しでjuliaを実行しよう。そうすると、下のようにRみたいなJuliaのREPLが立ち上がるよ!

$ julia
               _
   _       _ _(_)_     |  A fresh approach to technical computing
  (_)     | (_) (_)    |  Documentation: http://docs.julialang.org
   _ _   _| |_  __ _   |  Type "?help" for help.
  | | | | | | |/ _` |  |
  | | |_| | | | (_| |  |  Version 0.5.0 (2016-09-19 18:14 UTC)
 _/ |\__'_|_|_|\__'_|  |
|__/                   |  x86_64-apple-darwin14.5.0

julia> 1 + 2
3

julia>

ヘルプの調べ方は?

RみたいにREPLに?と打ってみよう。そうするとプロンプトがhelp?>に変わるから、知りたい関数の名前をそこに打ってEnterを押そう!

help?> sum
search: sum sum! sumabs summary sumabs2 sumabs! sum_kbn sumabs2! cumsum cumsum! consume cumsum_kbn isnumber isalnum SlotNumber

  sum(itr)

  Returns the sum of all elements in a collection.

  sum(A, dims)

  Sum elements of an array over the given dimensions.

  sum(f, itr)

  Sum the results of calling function f on each element of itr.

キーワード検索をする場合は、help?>プロンプトでキーワードを"..."で囲んで入力しよう!

help?> "variance"
Base.var
Base.cov
Base.varm

RStudioみたいなのはあるの?

RStudioほどじゃないかもしれないけど、JunoというAtomエディターをベースにした開発環境が人気だよ。他にはEclipse向けのプラグイン(https://github.com/JuliaComputing/JuliaDT)やVS Code向けのプラグイン(https://github.com/JuliaEditorSupport/julia-vscode)もあるよ。

EmacsVimの人は?

Emacsプラグイン(https://github.com/JuliaEditorSupport/julia-emacs)とVimプラグイン(https://github.com/JuliaEditorSupport/julia-vim)ももちろんあるよ。

ノートブックはあるの?

IJulia.jlを使って、Jupyter Notebookを利用できるよ。すぐ試してみるのにJuliaBoxを使うといいんじゃないかな?

CRANみたいなパッケージのレポジトリはあるの?

あるよ!パッケージの一覧がhttp://pkg.julialang.org/から閲覧できるよ。

パッケージのインストールはどうするの?

Rのinstall.packagesに当たるのがPkg.addだよ。もちろん依存解決もしてくれるよ。

GitHubにしかないパッケージのインストールは?

Pkg.cloneにレポジトリのURLを入れると野良パッケージでもインストールできるよ。ここでも依存解決をしてくれるよ。

Rを呼び出したいんだけど?

RCall.jlを使おう! RCall.jlはJuliaのコードにRを埋め込んだりRのREPLを使ったりできるよ。

Pythonも呼び出したいんだけど?

PyCall.jlを使おう! ついでに、C++が呼び出せるCxx.jlとかMATLABが呼び出せるMATLAB.jlなんてのもあるよ。

CやFortranのコードの呼び出しは?

JuliaはCやFortranのコードを呼び出す仕組みを標準で用意しているよ。 詳しくは後で説明するよ。

Rのサンプルデータを読みたいんだけど?

RDatasets.jlにRでよく使われるデータセットがパッケージされているよ!

ワークスペースを保存するにはどうするの?

JLD.jl@saveマクロを使おう。

julia> using JLD

julia> n = 5
5

julia> x = randn(n, n)
5×5 Array{Float64,2}:
 -2.18618     0.0831716   2.19386    0.806268  1.11444
 -0.99689    -0.187922   -0.138358  -0.141601  0.19058
  0.0971969   0.149858   -1.19826    1.15507   1.0969
 -0.140956   -0.727159    0.780267   1.17143   0.979918
 -0.161376   -0.162103   -0.175158  -0.402853  0.916248

julia> @save "workspace.jld"

ワークスペースを復元するには?

同じくJLD.jlの@loadマクロを使おう。

julia> using JLD

julia> @load "workspace.jld"
3-element Array{Symbol,1}:
 :ans
 :n
 :x

julia> n
5

julia> x
5×5 Array{Float64,2}:
 -2.18618     0.0831716   2.19386    0.806268  1.11444
 -0.99689    -0.187922   -0.138358  -0.141601  0.19058
  0.0971969   0.149858   -1.19826    1.15507   1.0969
 -0.140956   -0.727159    0.780267   1.17143   0.979918
 -0.161376   -0.162103   -0.175158  -0.402853  0.916248

コミュニティ

Juliaのコミュニティはどんな感じなの?

Juliaのコミュニティは小さいけれど、とてもオープンなコミュニティを形成しているよ。 特に初心者に優しい(Juliaユーザはだいたい初心者)ので、気軽にコミュニケーションを取れるよ。

質問はどこですればいいの?

英語ならDiscourseかStack Overflowだね。Stack Overflowで質問する時は"julia-lang"のタグをつけるようにしよう。 Discourseの方がJuliaハッカーが見てる率が高い気がするから、Juliaに限定した質問ならDiscourseの方がオススメかも。

日本語話者のコミュニティは?

日本ではJulia Tokyoが知る限り唯一のコミュニティかな? Facebook・Slack・GitHubのアカウントがあるし、不定期でミートアップも開催してるよ。

他にはどんなコミュニティがあるの?

Juliaのコミュニティは各ドメインに分かれてゆるくつながっているよ。 コミュニティのリストはここ(http://julialang.org/community/)にあるよ。 いくつか挙げると以下のドメインのコミュニティは結構活発かな。

ドメイン GitHub
統計 JuliaStats
最適化 JuliaOpt
データベース JuliaDB
GPU JuliaGPU
生物学 BioJulia
微分 JuliaDiff
グラフィクス JuliaGraphics
計量経済学 QuantEcon

これらで活発に開発されているパッケージは概ね安心して使えるかな。

どんな人達がJuliaを使ってるの?

Juliaの使い方は色々あるけれど、やっぱり数値計算目的に使っているユーザが多いかな。 Julia Computing, Inc.の事例を見ると、金融・経済分野を始めとしたデータ分析界隈で使われ始めている印象があるな。

HPCでの事例や科学論文でJuliaを使った計算の実装も最近よく見るようになってきたよ。

他には、StanfordやMITでは技術計算の授業にもJuliaが使われているよ。

トークとかの映像は?

YouTubeのチャンネルにカンファレンスの動画などがたくさんあるよ!

https://www.youtube.com/user/JuliaLanguage

Julia界のHadley Wickhamは?

Julia界には(良くも悪くも)あれくらい目立っている人はいないかな。 Juliaの重要なパッケージは各organizationがコミュニティとして開発・管理してるよ。

文法

Rの構文との対応関係を教えて?

いいよ!

代入

# R
x <- 100
# Julia
x = 100

分岐

# R
if (x == 0) {
    print("zero")
} else if (x < 0) {
    print("negative")
} else {
    print("positive")
}
# Julia
if x == 0
    println("zero")
elseif x < 0
    println("negative")
else
    println("positive")
end

for文

# R
for (i in 1:10) {
    print(i)
    if (i > 5) {
        break
    }
}
# Julia
for i in 1:10
    println(i)
    if i > 5
        break
    end
end

関数

# R
add <- function (x, y) x + y

add <- function (x, y) {
    return (x + y)
}
# Julia
add(x, y) = x + y

function add(x, y)
    return x + y
end

ライブラリ読み込み

# R
library(ggplot2)
# Julia
using DataFrames

たまにJuliaのコードにある@は何?

@から始まるコードはJuliaのマクロ呼出しだよ。 マクロはJuliaのコードを別のコードに書き換えたりするメタプログラミングの一種だよ。 例えば次のコードに出てくる@inboundsは配列アクセスの境界チェックを無くして少し高速化することができるし、@showはその時の値をいい感じに表示してくれるよ。

x = randn(10)
s = 0
@inbounds for i in 1:endof(x)
    s += x[i]
end
@show s

たまにJuliaのコードにあるr"..."は何?

それもJuliaのマクロだよ! 標準ライブラリではr"..."と書くことで文字列でなく正規表現を作ることができるよ。

julia> r"foo"
r"foo"

julia> typeof(r"foo")
Regex

julia> b"foo"
3-element Array{UInt8,1}:
 0x66
 0x6f
 0x6f

julia> typeof(b"foo")
Array{UInt8,1}

この仕組はユーザからも拡張可能なので、例えばBio.jlではDNA配列を作ったりするのに使ってるよ。

julia> using Bio.Seq

julia> dna"ACGTAG"
6nt DNA Sequence:
ACGTAG

データ

主な数値型は?

整数はInt型、倍精度浮動小数点数Float64型がデフォルトでよく使われるよ。 Intは32-bit環境なら32-bit、64-bit環境なら64-bitで表現されるよ。

Juliaは数値の型も豊富で、8-bitから128-bitまで符号の有り無しの組み合わせがすべて用意されているし、複素数(Complex{T})や有理数(Rational{T})もあるよ。

TRUEFALSEは?

Juliaでは全部小文字のtruefalseだよ。

欠損値は扱えるの?

JuliaではNullable{T}という型があって、型Tの欠損値を表すよ。 nothingという値もあるけどパフォーマンス上の理由であまりお勧めしないよ。

JuliaもRみたいに常に配列として数値を扱うの?

Juliaでは単一の値と配列は区別されるよ。例えば、3と書いたら整数の3を意味するし、[3]と書いたら1つの整数3からなるベクトルだよ。

3と書いたら浮動小数点数じゃなくて整数なの?

そうだよ、Rと違ってJuliaでは整数として扱われるよ。

じゃぁ浮動小数点数の3を作るには?

3.0と書けば浮動小数点数として扱われるよ。

どうやって確認するの?

Rみたいにtypeofという関数を使おう!

julia> typeof(3)
Int64

julia> typeof(3.0)
Float64

R> typeof(3)
[1] "double"

R> typeof(3.0)
[1] "double"

演算子はどんな感じ?

ほとんどRと同じだよ。

julia> 3 + 4
7

julia> 3 - 4
-1

julia> 3 * 4
12

julia> 3 / 4
0.75

R> 3 + 4
[1] 7

R> 3 - 4
[1] -1

R> 3 * 4
[1] 12

R> 3 / 4
[1] 0.75

Bool演算(&&, ||, !)も同じだよ。排他的論理和はv0.5では$だけど次期バージョンでxorという名前に変更される予定だよ。

listはあるの?

Juliaで一番近いのはDictかな? Pythondictと同じハッシュを使った連想配列だよ。

julia> d = Dict("one" => 1, "two" => 2, "three" => 3)
Dict{String,Int64} with 3 entries:
  "two"   => 2
  "one"   => 1
  "three" => 3

julia> d["two"]
2

data.frameは?

残念ながらJuliaの標準ライブラリにはないよ。でもDataFrames.jlというパッケージがJuliaのデータフレームの標準的な実装だよ。 この辺は後で詳しく述べるよ。

多次元行列も扱えるの?

もちろん扱えるよ! Rのmatrixみたいなことが次のようにできるよ。

julia> [1 2 3
        4 5 6]
2×3 Array{Int64,2}:
 1  2  3
 4  5  6

julia> [1 2 3; 4 5 6]  # セミコロンは改行と同じ扱い
2×3 Array{Int64,2}:
 1  2  3
 4  5  6

R> matrix(c(1, 2, 3, 4, 5, 6), nrow=2)
     [,1] [,2] [,3]
[1,]    1    3    5
[2,]    2    4    6

空の配列の作成は?

ベクトル(1次元配列)なら[]だよ。

julia> []
0-element Array{Any,1}

ただ要素の型が分かっている場合はそれを教えてあげるとJuliaは高速に動くよ。

julia> Int[]
0-element Array{Int64,1}

julia> Float64[]
0-element Array{Float64,1}

matrix(0, 3, 4)みたいに0で初期化された配列の作成は?

zeros関数を使おう!

julia> zeros(3, 4)  # 3行4列の行列
3×4 Array{Float64,2}:
 0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0

julia> zeros(Int, 3, 4)
3×4 Array{Int64,2}:
 0  0  0  0
 0  0  0  0
 0  0  0  0

1で初期化された配列ならonesだよ。

配列要素へのアクセスはどうするの?

基本的にRと同じだよ。Juliaの配列も1始まりでcolumn majorだよ。 ただ、Rだとx[i,]とかx[,j]と書くところをJuliaではx[i,:]とかx[:,j]とか書くよ。

julia> x = [1 2 3; 4 5 6]
2×3 Array{Int64,2}:
 1  2  3
 4  5  6

julia> x[1,2]
2

julia> x[2,2:3]
2-element Array{Int64,1}:
 5
 6

julia> x[:,2]
2-element Array{Int64,1}:
 2
 5

2:3みたいなのはベクトル?

Juliaではm:nmからnまでの範囲(Range)を表すよ。

julia> 2:3
2:3

julia> typeof(2:3)
UnitRange{Int64}

ベクトルと違ってメモリーを消費しないので、長い範囲でも一瞬で作れるよ。

characterみたいな文字列型は?

JuliaではString型があるよ。

文字列はUnicodeも使える?

使えるよ! JuliaではUTF-8エンコーディングしてあるよ。

文字列操作関数を教えて?

stringrパッケージの関数と大雑把に対応を付けたよ。

R (stringr) Julia
str_c string
str_detect ismatch
str_extract match
str_dup ^
str_split split
str_replace replace
str_trim strip, lstrip, rstrip
str_pad lpad, rpad
str_to_lower lowercase
str_to_upper uppercase

as.doubleとかas.integerみたいな型変換はどうするの?

Juliaの型変換は基本的にconvert関数を使うよ。

julia> convert(Float64, 42)
42.0

julia> convert(Int64, 42.0)
42

ただし文字列をパースして数値にするにはparseを使うよ。

julia> parse(Int, "42")
42

julia> parse(Float64, "42")
42.0

行名や列名はどうつけるの?

Juliaの標準では行名や列名をサポートしていないよ。 でもNamedArrays.jlという配列に名前を付けられるようにするパッケージはあるよ。

欠損値のある配列はどう扱えばいいの?

欠損値の扱いはまだJulia界隈でも完全な合意に達してない課題なので、ちょっとむずかしいな。 前にも述べたようにNullable{T}というデータ型があるけど、これに関する操作が将来変わる可能性があるよ。 でも、基本的にはNullableArrays.jlを使おう!

因子(factor)を扱うにはどうするの?

これもまだJuliaではちょっと扱いが難しいけど、CategoricalArrays.jlというパッケージがあるよ。 NullableArrays.jlやCategoricalArrays.jlは近いうちにDataFrames.jlに採用される予定のようだよ。

関数

関数はどうやって定義するの?

functionキーワードを使ったスタイルと短縮記法があるよ。

# 通常の記法
function add(x, y)
    return x + y
end

# 短縮記法
add(x, y) = x + y

デフォルト引数やキーワード引数は?

あるよ! でもJuliaではデフォルト引数とキーワード引数で記法が分かれているよ。

# デフォルト引数
function default(x, y=1, z=10)
    return x + y + z
end

# キーワード引数 (セミコロンに注目)
function keyword(x; y=1, z=10)
    return x + y + z
end

試しに実行してみよう

julia> default(1)
12

julia> default(1, 0)
11

julia> default(1, 0, 5)
6

julia> keyword(1)
12

julia> keyword(1, z=5)
7

julia> keyword(1, z=5, y=0)
6

クロージャも作れる?

もちろん!

function counter()
    i = 0
    function x(n=1)
        i += n
        return i
    end
end
julia> c = counter()
(::x) (generic function with 2 methods)

julia> c()
1

julia> c()
2

julia> c()
3

julia> c(3)
6

julia> c()
7

+とかは関数じゃないの?

関数だよ! 前置記法で書くこともできるよ!

julia> +(1, 2)
3

R> `+`(1, 2)
[1] 3

じゃぁRみたいにifとかforも関数なの?

Juliaではifforwhileは関数じゃないよ。

遅延評価はあるの?

Juliaには遅延評価はなくてすべて正格評価だよ。つまり関数に渡された引数は呼び出し前にかならず評価されるよ。

applyみたいな関数を配列の要素に適用する関数は?

mapを使おう! ただし関数が最初の引数だよ。

julia> map(log, [1,2,3])
3-element Array{Float64,1}:
 0.0
 0.693147
 1.09861

julia> map(exp, [1,2,3])
3-element Array{Float64,1}:
  2.71828
  7.38906
 20.0855

R> sapply(c(1,2,3), log)
[1] 0.0000000 0.6931472 1.0986123

R> sapply(c(1,2,3), exp)
[1]  2.718282  7.389056 20.085537

ただ最近はブロードキャストを使う方法が一般的なので、こちらを使ったほうが良いよ。

julia> log.([1,2,3])
3-element Array{Float64,1}:
 0.0
 0.693147
 1.09861

julia> exp.([1,2,3])
3-element Array{Float64,1}:
  2.71828
  7.38906
 20.0855

この関数名の後のドットは何?

ブロードキャスト関数呼び出しの構文糖衣だよ。例えばf.(x)broadcast(x, f)に変換されるよ。 broadcastmapの凄い版だよ。

ベクトル同士の足し引きとかは?

+-がそのまま使えるよ。でも.を付けた.+.-なら左右で配列のサイズが違った場合でもRのように動いてくれるよ。

julia> [1,2,3] + [4,5,6]
3-element Array{Int64,1}:
 5
 7
 9

julia> [1,2,3] .+ [4]
3-element Array{Int64,1}:
 5
 6
 7

julia> [1,2,3] .+ 4
3-element Array{Int64,1}:
 5
 6
 7

R> c(1, 2, 3) + c(4, 5, 6)
[1] 5 7 9

R> c(1, 2, 3) + c(4)
[1] 5 6 7

R> c(1, 2, 3) + 4
[1] 5 6 7

書いた関数が思ったより遅いんだけど?

一番よくある原因が関数内で使われている変数の型が不安定になっているせいだよ。 例えば次のコードはxsの要素がIntのときには高速に動くけど、浮動小数点数の配列とかを渡すとすごく遅いよ。

function sumup(xs)
    s = 0
    for x in xs
        s += x
    end
    return s
end
julia> xs = rand(1:10, 10000);

julia> sumup(xs);  # ウォームアップ

julia> @time sumup(xs)
  0.000011 seconds (5 allocations: 176 bytes)
54649

julia> xs = rand(10000);

julia> sumup(xs);  # ウォームアップ

julia> @time sumup(xs)
  0.000337 seconds (30.00 k allocations: 468.906 KB)
5015.030044010504

原因は変数sは整数型で初期化されてるのに、ループの中で浮動小数点数と足すことでsの型が整数か浮動小数点数かJuliaには分からなくなってしまうせいだよ。これをJulia界隈では型の不安定さ(type instability)と言ってるよ。 これを回避するにはeltypezero関数を次のように使うことで、sの型をxsに合わせて設定できるよ。

function sumup(xs)
    s = zero(eltype(xs))
    for x in xs
        s += x
    end
    return s
end

こうすることでさっきより何十倍も速くなったね!

julia> sumup(xs);

julia> @time sumup(xs)
  0.000015 seconds (5 allocations: 176 bytes)
5015.030044010504

関数のプロファイルを取るにはどうするの?

@profileマクロと、Profile.print()をあわせて使おう! @profile <code>でコード実行のプロファイルをとって、Profile.print()でプロファイル結果をプリントできるよ。

C言語の関数の呼び出し方法を教えて!

JuliaにはccallというC言語を呼び出す専用の命令があるよ。例えば、GNU libcに実装されているexp関数を呼び出すなら次のように書けるよ。

julia> ccall((:exp, "libc"), Float64, (Float64,), 1.0)
2.718281828459045

julia> ccall((:exp, "libc"), Float64, (Float64,), -1.0)
0.36787944117144233

ccallは最初に呼び出す関数を指定して、その後に返り値の型と引数の型を指定して、最後に引数を渡すよ。上の例で言えば、(:exp, "libc")が呼び出す関数で、Float64が返り値の型、(Float64,)が引数の型で、1.0-1.0がCの関数に渡される実際の引数だよ。

Juliaは数値以外にも文字列やより複雑な構造体をCのライブラリとやり取りすることができるように設計されているので、より詳しくはマニュアルのここを読んでみてね。 具体的な例としてはJuliaの標準ライブラリやLibz.jlなんかが勉強になるよ。

型システムとメソッド

Juliaの型システムはどうなってるの?

Juliaの型は大きく具体型(concrete type)と抽象型(abstract type)に分けられていて、メソッドの選択は多重ディスパッチというシステムを中心に構築されているよ。

具体型とか抽象型って何?

具体型はインスタンス化(オブジェクトを作れる)できる型で、抽象型はそうでない型だよ。 例えばFloat641.0など値を作れるけどRealという抽象型は値を作れないよ。

値を作れない抽象型は何の役に立つの?

複数の型をまとめるのに使われるよ。KeitaNakamuraさんのコードを使わせてもらうと、例えばRealを頂点とする型は次のような構造になっているよ。

julia> print_tree(Real)
Real
├─── AbstractFloat
|    ├─── BigFloat
|    ├─── Float16
|    ├─── Float32
|    └─── Float64
├─── Integer
|    ├─── BigInt
|    ├─── Bool
|    ├─── Signed
|    |    ├─── Int128
|    |    ├─── Int16
|    |    ├─── Int32
|    |    ├─── Int64
|    |    └─── Int8
|    └─── Unsigned
|         ├─── UInt128
|         ├─── UInt16
|         ├─── UInt32
|         ├─── UInt64
|         └─── UInt8
├─── Irrational{sym}
└─── Rational{T<:Integer}
function print_tree(typ)
    println(typ)
    print_tree(typ, [])
end

function print_tree(typ, level)
    stypes = subtypes(typ)
    for stype in stypes
        if stype !== stypes[end]
            println(join(level) * "├─── $stype")
            push!(level, "|    ")
        else
            println(join(level) * "└─── $stype")
            push!(level, "     ")
        end
        print_tree(stype, level)
        pop!(level)
    end
end

クラス(型)の定義はどうやるの?

RのS4に似た仕組みがJuliaにはあるよ。例えば以下の2つの定義は大体対応しているよ。

# R
setClass("Person", representation(name = "character", age = "numeric"))
# Julia
type Person
    name::String
    age::Int
end

Juliaはデフォルトコンストラクタを準備しているから、オブジェクトは次のように作成できるよ。

julia> hadley = Person("Hadley", 31)
Person("Hadley",31)

julia> hadley.name
"Hadley"

julia> hadley.age
31

ただし、Juliaのtypeで宣言したオブジェクトは変更可能なのでその振る舞いはRの参照クラス(RC)に近いかな。

継承は?

JuliaはsetClasscontains引数に当たるような継承の仕組みはないよ。つまり、複数の型間で構造を継承する方法はないよ。

その代わり、定義する型をある抽象型のサブタイプにする仕組みはあって、Juliaではこちらを使うよ。 例えば、Juliaで要素の数が3つのベクトルっぽい型を作る場合には、次のようにAbstractVectorのサブタイプとして宣言すると良いよ。

type Triplet{T} <: AbstractVector{T}
    x::T
    y::T
    z::T
end

ここで、Tと書かれているのは型パラメータで、T = Int64T = Float64が代入されて具体的な型になるよ。 また、各抽象型には実装していることが期待されるメソッドがあるので、それを実装するべきだよ (AbstractVectorならBase.sizeBase.getindexなどを実装すべき)。

メソッドの定義は?

通常の関数定義に型指定したものがRのsetMethodに近いよ。

# R
setGeneric("greet", function(x) standardGeneric("greet"))
setMethod("greet", signature("Person"), function(x) paste("Hello, ", x@name))
# Julia
function greet(x::Person)
    return string("Hello, ", x.name)
end

型指定は複数に引数に対しても指定できるし、抽象型も指定できるよ。これをJuliaでは多重ディスパッチと言うよ。

前処理

CSV/TSVファイルの読み込みは?

3通りの方法があるよ!

  1. 標準ライブラリのreadcsv
  2. DataFrames.jlのreadtable
  3. CSV.jlCSV.read

どれがいいの?

標準ライブラリのreadcsvは標準なのですぐ使えるけど、行列としてデータを読み込むから、場合によってはちょっと使いづらいかな? DataFrames.jlのreadtableはデータフレームを読み込むけど、将来的に使われなくなる可能性が高いよ。 なので、現在ではCSV.jlのCSV.readを使うのがオススメだよ。

julia> using DataFrames

julia> using CSV

julia> head(CSV.read("iris.csv"))
6×5 DataFrames.DataFrame
│ Row │ SepalLength │ SepalWidth │ PetalLength │ PetalWidth │ Species  │
├─────┼─────────────┼────────────┼─────────────┼────────────┼──────────┤
│ 1   │ 5.1         │ 3.5        │ 1.4         │ 0.2        │ "setosa" │
│ 2   │ 4.9         │ 3.0        │ 1.4         │ 0.2        │ "setosa" │
│ 3   │ 4.7         │ 3.2        │ 1.3         │ 0.2        │ "setosa" │
│ 4   │ 4.6         │ 3.1        │ 1.5         │ 0.2        │ "setosa" │
│ 5   │ 5.0         │ 3.6        │ 1.4         │ 0.2        │ "setosa" │
│ 6   │ 5.4         │ 3.9        │ 1.7         │ 0.4        │ "setosa" │

詳しい使い方はCSV.jlのヘルプを参照してね。

SQLiteのデータを読みたいんだけど?

SQLite.jlを使おう。SQLを発行した結果はDataFrameとして返してくれるよ。

JSONを扱うには?

JSON.jlを使おう。

julia> using JSON

julia> JSON.parse("""{"foo": 100, "bar": [1.1, 2.0]}""")
Dict{String,Any} with 2 entries:
  "bar" => Any[1.1,2.0]
  "foo" => 100

XMLやHTMLのデータを処理するには?

拙作のEzXML.jlがオススメだよ。長く使われているLightXML.jlというのもあるけれど、機能が少ないし、あまりインターフェースがきれいじゃないね...

julia> using EzXML

julia> readxml("ex1.xml")
EzXML.Document(EzXML.Node(<DOCUMENT_NODE@0x00007fcda9030400>))

julia> readxml("ex1.xml") |> print
<?xml version="1.0" encoding="UTF-8"?>
<bookstore>
  <book category="COOKING" tag="first">
    <title lang="en">Everyday Italian</title>
    <author>Giada De Laurentiis</author>
    <year>2005</year>
    <price>30.00</price>
  </book>
  <book category="CHILDREN">
    <title lang="en">Harry Potter</title>
    <author>J K. Rowling</author>
    <year>2005</year>
    <price>29.99</price>
  </book>
</bookstore>

でもHTMLに関してはGoogleのgumboをラップしたGumbo.jlの方が良いかな?

dplyrは?

残念ながらdplyrほどの完成度のパッケージは今のところJuliaには無いかな。 ただJuliaの人にはR使いも多いので、dplyrの利便性は誰もが認めるところだね。 それで最近は似たことをJuliaでやろうとしている人達が出てきて、いくつかパッケージができているよ。

Juliaのデータフレームの近況と将来についてはここにまとまってるよ。

UTF-8以外のファイルを読み込みたいんだけど?

StringEncodings.jlを使おう! libiconvを使ってるのでSJISみたいなアレなエンコーディングでもちゃんと扱えるよ!

統計/機械学習

平均とか分散の計算は?

mean, median, var, cov, corはRと同名の関数が標準ライブラリにあるよ。標準偏差sdでなくstdになってるので注意だよ。

rowSumsとかcolSumsとかは?

sum(array, dim)が使えるよ。

julia> x
2×3 Array{Int64,2}:
 1  2  3
 4  5  6

julia> sum(x, 1)  # colSums
1×3 Array{Int64,2}:
 5  7  9

julia> sum(x, 2)  # rowSums
2×1 Array{Int64,2}:
  6
 15

他にもmeanvarなんかでも同じことができるよ。

乱数生成はどうするの?

一様乱数runif(n)に対応するのはrand(n)正規分布rnorm(n)に対応するのはrandn(n)だよ。

ガンマ分布とかポアソン分布みたいな他の分布は?

Distributions.jlを使おう。

PCA(主成分分析)をしたいんだけど?

MultivariateStats.jlを使おう。

julia> using MultivariateStats

julia> using RDatasets

julia> iris = dataset("datasets", "iris")

julia> pca = fit(PCA, Matrix(iris[:,1:4])')  # 行列に変換して転置
PCA(indim = 4, outdim = 3, principalratio = 0.99479)

julia> transform(pca, Matrix(iris[:,1:4])')
3×150 Array{Float64,2}:
 2.68413     2.71414    2.88899     2.74534     2.72872     2.28086   …  -1.94411   -1.52717   -1.76435    -1.90094   -1.39019
 0.319397   -0.177001  -0.144949   -0.318299    0.326755    0.74133       0.187532  -0.375317   0.0788589   0.116628  -0.282661
 0.0279148   0.210464  -0.0179003  -0.0315594  -0.0900792  -0.168678     -0.177825   0.121898  -0.130482   -0.723252  -0.36291

GML(一般化線形モデル)を扱いたいんだけど?

GLM.jlを使おう。Rでも有名なDouglas Bates先生の作ったパッケージだから安心だね!

Lassoは?

Julia実装のLasso.jlかglmnetをラップしたGLMNet.jlが良いと思うよ。

Stanを使いたいんだけど?

Stan.jlを使おう!

行列積の計算は?

Rの%*%がJuliaでは*だよ。

julia> A = randn(2, 3);

julia> B = randn(3, 4);

julia> A * B
2×4 Array{Float64,2}:
 2.42828   -1.39917   0.28215   -1.61981
 0.292669  -0.411146  0.234041   1.82686

R> A <- matrix(rnorm(6), 2)

R> B <- matrix(rnorm(12), 3)

R> A %*% B
          [,1]       [,2]       [,3]      [,4]
[1,] 0.8448714 -0.6084743 -0.7458281 -0.653775
[2,] 0.2133118 -0.9305351  0.3933490  1.105863

solveみたいな線形方程式の解の計算は?

A * X = Bを解くには、A \ Bと計算するよ。

julia> A = randn(3, 3)
3×3 Array{Float64,2}:
  0.864382  -0.945314  -0.496754
  0.119709   0.83993   -0.962021
 -0.737123   0.221293   0.205341

julia> X = randn(3, 2)
3×2 Array{Float64,2}:
 -0.224859   1.64374
  0.0726233  0.898944
 -0.387804   0.697085

julia> B = A * X
3×2 Array{Float64,2}:
 -0.0703728   0.224752
  0.407157    0.28121
  0.102188   -0.869566

julia> A \ B
3×2 Array{Float64,2}:
 -0.224859   1.64374
  0.0726233  0.898944
 -0.387804   0.697085

julia> A \ B ≈ X
true

Aが実対称行列であることが分かってるんだけど?

cholfact(A)でコレスキー分解してからでも同様に計算できるよ!

julia> A = randn(3, 3);

julia> A = A'A
3×3 Array{Float64,2}:
  1.13074     2.65079   -0.0602727
  2.65079     7.18991   -0.474305
 -0.0602727  -0.474305   0.488933

julia> issymmetric(A)
true

julia> X = randn(3, 2)
3×2 Array{Float64,2}:
  0.809623   0.595824
 -0.25161   -2.33654
  0.97426    1.46868

julia> B = A * X;

julia> A \ B
3×2 Array{Float64,2}:
  0.809623   0.595824
 -0.25161   -2.33654
  0.97426    1.46868

julia> cholfact(A) \ B
3×2 Array{Float64,2}:
  0.809623   0.595824
 -0.25161   -2.33654
  0.97426    1.46868

log(sum(exp(x)))とかlog(exp(x)+1)とかを安定的に計算したいんだけど?

StatsFuns.jlに標準ライブラリにはないけど、数値計算でよく見る関数の実装があるからそっちを探してみよう。

julia> x = randn(10)
10-element Array{Float64,1}:
  0.0410153
 -0.407537
  1.26758
 -0.52522
 -0.0969684
 -2.52716
 -0.633663
 -1.2695
 -0.474312
 -0.800294

julia> log(sum(exp(x)))
2.165782421890441

julia> logsumexp(x)
2.165782421890441

julia> log(sum(exp(1000x)))  # 1000倍すると計算できない
Inf

julia> logsumexp(1000x)      # StatsFuns.jlの実装なら大丈夫
1267.5802777924926

プロット作るにはどうするの?

残念ながら、今のところJuliaには決定的なプロットのための仕組みがあるわけではないので、次のものを参考に好みに合ったのを使おう!

  • Gadfly.jl
    • 最も歴史のあるJuliaの
    • 一番スターが多い
  • Plots.jl
    • 新興で急成長中のライブラリ
    • 開発が活発
  • PlotlyJS.jl
  • PyPlot.jl
    • Pythonのmatplotlibのラッパーライブラリ
    • 最も高機能かも
  • GR.jl
  • UnicodePlots.jl
    • Unicode文字を使ってターミナル上でプロットできるライブラリ
    • ごく簡単なプロットをチラッと見たいときに向いてる

関数の最大・最小値を探したいんだけど?

Optim.jlを使おう!

julia> using Optim

julia> rosenbrock(x) = (1.0 - x[1])^2 + 100.0 * (x[2] - x[1]^2)^2
rosenbrock (generic function with 1 method)

julia> optimize(rosenbrock, randn(2))
Results of Optimization Algorithm
 * Algorithm: Nelder-Mead
 * Starting Point: [0.4289709134439212,-1.3041568238701216]
 * Minimizer: [1.0000314591723356,1.0000648302370259]
 * Minimum: 1.354834e-09
 * Iterations: 69
 * Convergence: true
   *  √(Σ(yᵢ-ȳ)²)/n < NaN: false
   * Reached Maximum Number of Iterations: false
 * Objective Function Calls: 136

R> rosenbrock <- function (x) (1 - x[1])^2 + 100 * (x[2] - x[1]^2)^2

R> optim(rnorm(2), rosenbrock)
$par
[1] 0.9995882 0.9991591

$value
[1] 2.002849e-07

$counts
function gradient
     129       NA

$convergence
[1] 0

$message
NULL

Juliaの機械学習パッケージはどうなってるの?

Juliaユーザには機械学習を常用する人も多いけれど、PythonやRほど充実はしてないよ。 一応、ScikitLearn.jlというscikit-learnのラッパーがJuliaのパッケージになってるので、scikit-learnでできることはできるはずだよ。

Gradient boostingとかは?

XGBoost.jlが使えるよ! 名前から分かる通りXGBoostへのJuliaのインターフェースだよ。

Deep learningのライブラリは?

MXNet.jlというMXNetのパッケージがあるよ。またDMLCだね。他にはMocha.jlというJulia製のフレームワークもあるけれど、メインの開発者がmxnetの方に移ったようなので今はそれほど開発が活発でもないかな。 他にはTensorFlow.jlというPythonのTensorFlowのラッパーライブラリもあるよ。

もっとJuliaっぽいやつだと、Transformations.jlとかNAISTの進藤先生が作っているMerlin.jlというパッケージもあるよ。

結局、Juliaでデータ解析はどうなの?

Juliaは機械学習アルゴリズムを組むには良い言語だと思うけど、ライブラリやツールもとても少ないね。 なので現状、インタラクティブな探索的データ解析(exploratory data analysis)をするような場合にはJuliaよりRの方をオススメするよ。

ただ、最近JuliaMLという機械学習に特化したコミュニティができてJuliaで気軽に機械学習ができるようにしていたり、データフレームに本格的に見直そうとしていたりとか、改善に向けて動いているよ。今ならやることがたくさんあるので、活動に参加してみてはどうかな?

Julia現状確認 (言語編)

julia

これはJuliaアドベントカレンダー1日目の記事です。2016年におけるJuliaの「現状」を確認していこうと思います。

Julia言語とは

Juliaとは一体どのような言語でしょうか。julialang.orgの最初の説明を用いると、「他の技術計算環境のユーザに馴染みのある構文を備えた、高レベル・高パフォーマンスな技術計算にための動的プログラミング言語」("a high-level, high-performance dynamic programming language for technical computing, with syntax that is familiar to users of other technical computing environments.")です。ここで言う「他の技術計算環境」とは、MATLABPythonのことを指していると考えられます。実際、Juliaの構文や関数名はMATLABとの類似性が高いため、Octaveのようなオープンソース版のMATLAB代替と捉えられることも多いようですが、プログラミング言語としてはかなり性質が異なります。また、Pythonとは文法や関数名の面ではそれほど類似性は無いように思いますが、似ている点も多数あります。

Juliaの最大の特徴を挙げるならば、動的プログラミング言語でありながらC言語などに匹敵するほど実行が高速であるという点でしょう。多くの技術計算では関数の呼び出しやループが頻繁に用いられますが、Juliaでは実行時コンパイル(JIT)のおかげでこれらの実行コストがネイティブコードにコンパイルするプログラミング言語並になっています。しかし、Juliaが行う実行時のコンパイルは、他の言語でよくみられるJITコンパイラとは性質を異にするようです。PyPy(PythonJITコンパイラ実装)やLuaJIT(LuaJITコンパイラ実装)はTracing JITと呼ばれるもので、実行時にプロファイルを行ってホットスポットを同定し、実行時の情報を集め、コンパイルと最適化を行います。しかし、JuliaのJITコンパイラはコードの実行前に型推論と最適化を行い、その後コードを実行します。型推論の結果、型が決定できる部分では強力な最適化を行い、曖昧性が残る部分では実行時の型に合わせて実行できるようにしています。

Julia言語の現状

Julia言語の実装であるJulia処理系はすべてGitHub上で開発されています。レポジトリはこちら( https://github.com/JuliaLang/julia)です。2016年12月1日現在の最新のリリース版はv0.5.0で、9月20日頃リリースされたものです。次期バージョンはv0.6でほぼ間違いないと思います。気になるのはv1.0のリリース時期ですが、今年のJuliaConでJuliaの生みの親のひとりであるStefan Karpinski氏が宣言したことによると、Julia 1.0のリリースは2017年を予定しているようです。次期バージョンのv0.6はv1.0のα版的な意味合いで開発されているようなので、v0.6が出てすぐv1.0がリリースされるものと思われます。しかし実は、v1.0が出てしばらく1.0のままということはなく、1-2年のうちにv2.0も出すつもりのようです。

振り返って、現行のv0.5ではどのような変更があったのかを確認してみます。まず性能面で大きな影響があったのが、クロージャが高速になったことです。v0.4までのJuliaでは、関数がすべて一つの型にまとめられており最適化の恩恵を受けにくかったのですが、v0.5からは関数が別々の型を持つようになりました。これによって、無名関数などの関数の呼び出しが高速化され、安心して使えるようになりました。例えば、map(x -> 2x, array)のようなコードがforループと同等の速度が出るようになります。また、ASCIIとUTF-8に分かれていた文字列型がString型に統一され、内部のエンコーディングを気にせず使えるようになりました。これで、Julia内では常にUTF-8エンコーディングのみを考慮すればよく、プログラムを書くのが大変楽になりました。さらに、新機能として、実験的にマルチスレッディングのサポートが入りました。今のところ、Julia内部でまだマルチスレッディングに対応していない部分(例えばIOなど)があるので、常に問題なく使えるわけではありませんが、競合などがない場合にはforループに@threadsをつけるだけで中身を並列に実行してくれます。これらの変更の他にも、generator式(i for i in 1:10などが式として受け渡しできる)やfusing構文(sin.(cos.(x))などがベクトルxに対して一時配列を作らない)などの新しい機能が導入されています。

では、v0.6はどうなるでしょうか。v0.6はv1.0のリリースのためにいつもより早めにfeature freezeが入るようですので、そろそろ新機能や変更点などが出揃うと思います。その中で、おそらく一番大きな変更が型システムの変更になるでしょう。実装はここ(https://github.com/JuliaLang/julia/pull/18457)にプルリクがあるので、すぐに試すことができます。これを少し詳しく説明してみます。

JuliaにはSetと呼ばれる集合を表す型があります。これは型パラメータを取ることができる型で型パラメータをTとすると、Set{T}と書けます。つまりSet{T}T型の要素を持つ集合です。このとき型パラメータTはどのような型でも良いのですが(注: 実行時にはhashを実装しているなどの制約がある)、逆にTに何らかの制約(T<:Integerなど)を付けた型は作ることができませんでした。しかし、次期バージョンのJuliaの型システムでは、Set{T} where T <: IntegerなどのようにT型に制約をつけることができます。以下の例では、型制約がない場合とある場合でどのようにサブタイピングの結果が異なるかを示しています。

julia> Set{T} where T  # 型制約なし
Set{T} where T

julia> Set{Int} <: (Set{T} where T)
true

julia> Set{Float64} <: (Set{T} where T)
true

julia> Set{T} where T <: Integer  # 型制約あり
Set{T} where T<:Integer

julia> Int <: Integer
true

julia> Set{Int} <: (Set{T} where T <: Integer)
true

julia> Float64 <: Integer
false

julia> Set{Float64} <: (Set{T} where T <: Integer)
false

サブタイピングに使えるということは、多重ディスパッチにも使えるということです。実際、以下の様な引数の型も書けます。

julia> f(xs::(Tuple{V,V} where V <: AbstractVector{T}), x::T) where T <: Integer = (xs[1] + x, xs[2] - x)
f (generic function with 1 method)

julia> f(([1,2,3], [1,2,3]), 2)
([3,4,5],[-1,0,1])

これが実際どのように役立つのかは、実例が無いので私にもまだよく分かっていません。また、記法に関してはまだ変更が入るかもしれません。

もうひとつ重要な新機能は、依存する関数の自動更新です。Julia v0.5では以下のように関数gを実行した後fを更新してもfの変更は反映されませんでした。

julia> f(x) = x + 1
f (generic function with 1 method)

julia> g(x) = f(2x)
g (generic function with 1 method)

julia> g(1)
3

julia> f(x) = x - 1
WARNING: Method definition f(Any) in module Main at REPL[1]:1 overwritten at REPL[4]:1.
f (generic function with 1 method)

julia> g(1)
3

これでは、開発中など頻繁に関数を書き換えるときにいちいち再起動する必要があります。Juliaはコードをコンパイルするためこのような動的な変更は実装不可能だと私は思っていましたが、現在進行中の変更(https://github.com/JuliaLang/julia/pull/17057)を取り入れると、fの変更が正しくgの実行に反映されるようになります。

julia> f(x) = x + 1
f (generic function with 1 method)

julia> g(x) = f(2x)
g (generic function with 1 method)

julia> g(1)
3

julia> f(x) = x - 1
WARNING: Method definition f(Any) in module Main at REPL[1]:1 overwritten at REPL[4]:1.
f (generic function with 2 methods)

julia> g(1)
1

大変不思議です。どうやらworldという利用可能なメソッドをチェックできるようにする機構が組み込まれているようなのですが、私には詳しく分かりません。もし何か分かりましたら、また記事にしようと思います。

この記事ではJulia言語の現状を簡単にまとめてみました。Julia言語以外のコミュニティの動向などは、他の機会にまたご紹介しようと思います。明日は、Ishitaさんの担当です。

勉強会が当日キャンセルされた話

注意
  • 2016年10月27日 15時 追記
  • 2016年11月2日 16時 タイトル変更

TECH_SALONという団体が主催しているfreee株式会社のJulia言語の勉強会に参加しようとしたら、当日の午後に突如勉強会の開催が取り消されました。ここがconnpassのイベントページです。

techsalon.connpass.com

Julia言語の勉強会ということで、どんな人が参加するのか楽しみにしていただけに大変残念です。 イベントページには、参加者に向けてわざわざ次のような注意事項が記載されているにも関わらず、主催者側が当日にキャンセルするというのは理解しがたいものがあります。

やむを得ない事情でキャンセルされる方は必ず3日前までにキャンセル申請をしてください。それ以降にキャンセルされる方は今後のTECH_SALONの勉強会の参加をお断りする可能性があります。

主催者側からどのような理由があって中止になったのかは明らかにされていませんが、以下の様なメッセージが届きました。

会場提供をしてくださっているfreee株式会社様のご都合により、急遽開催を見送らせて頂く事となりました。直前の連絡で大変申し訳ありません。次回、開催の際、ご連絡致しますので今後ともに何卒宜しくお願い致します。


(追記 2016年10月27日 15時)

主催団体より、中止事由の引用部分を次の内容に修正して欲しいとの連絡がありました。

勉強会の主催を務めるTECH_SALON側の判断により、

実際にはfreee株式会社の判断ではなく、TECH_SALONという主催団体側の判断とのことです。


会場となるfreee株式会社の都合によりキャンセルされたようですが、何があったのでしょうか。私は、連絡があったのを知らず直接会社まで行ったのですが、特に変わった様子は見受けられませんでした。

勉強会の講師の人に問い合わせると、以下の様な返信が届きました。

参加人数があまりにも少なかったため急遽開催を取りやめたようです。しかし、これは既に数日前に分かっていたことではないのでしょうか。

きっと、私には理解できない何か深遠な理由が裏に隠されているのだと思いますが、TECH_SALONの勉強会はこのような主催自体のドタキャンが発生する可能性があるため、参加の際にはご注意下さい。

Julia言語の0.5の変更点

9月20日にJulia言語の最新版である0.5がリリースされました。Juliaのメーリングリストに投稿されたアナウンスメントはこちらです: https://groups.google.com/d/msg/julia-users/J2DiH1GnM8o/aO2Ku8o-CgAJ

きっと近いうちに本家のブログで詳しい変更点の紹介があると思いますが、私のブログでも一足先に主要な変更点をご紹介しようと思います。

クロージャの効率化

とりわけ重要な変更点として挙げられるのがJuliaのクロージャが効率化されたことです。 0.4までのJuliaでは、Juliaの関数に関数を渡したりJuliaで関数を返すような関数を作ると、その実効速度が極めて遅いことが問題でした。 これは、すべての関数が Functionという型にまとめられていたせいで、Juliaのコンパイラが特化したコードを吐けないせいでした。 これが0.5では関数は各々自分専用の型を持つように変更されたため、クロージャを使ったコードでも遅くなるということが無くなりました。

関数を受け取る関数の例としてよく上げられるのがmap関数でしょう。 以下の簡単なベンチマークを見てもJulia 0.5では0.4に比べて10倍ほど速度が向上しています。

Julia 0.5.0:

julia> f = x -> 2x
(::#1) (generic function with 1 method)

julia> x = rand(100000);

julia> map(f, x);

julia> @time map(f, x);
  0.000572 seconds (134 allocations: 789.078 KB)

Julia 0.4.6:

julia> f = x -> 2x
(anonymous function)

julia> x = rand(100000);

julia> map(f, x);

jjulia> @time map(f, x);
  0.005628 seconds (200.01 k allocations: 3.815 MB)

クロージャを使うようなプログラムではこれでデザインの幅が広がり、今まで実行効率の観点からできなかったコードが書けるようになります。

ジェネレータ式

ジェネレータ式という新しい式が0.5から加わりました。 これは、イテレータをラップして式にしたようなもので、値を次々生成することができる式です。 簡単な例では、 2x for x in 1:9 のようなものが普通の値として扱えるようになりました。

julia> g = (2x for x in 1:9)
Base.Generator{UnitRange{Int64},##5#6}(#5,1:9)

julia> collect(g)
9-element Array{Int64,1}:
  2
  4
  6
  8
 10
 12
 14
 16
 18

もちろん関数はこのジェネレータを値として受け取ることもできるため、次のような計算もできます。

julia> sum(g)
90

また、ifで値をフィルターすることも可能です。

julia> h = (2x for x in 1:9 if isodd(x))
Base.Generator{Filter{Base.#isodd,UnitRange{Int64}},##7#8}(#7,Filter{Base.#isod
d,UnitRange{Int64}}(isodd,1:9))

julia> collect(h)
5-element Array{Int64,1}:
  2
  6
 10
 14
 18

julia> sum(h)
50

文字列型の統合

Julia 0.4の標準ライブラリに大量にあった文字列型が整理され、0.5では主にString型とSubString型のみになりました。 String型は従来のUTF8String型に当たるもので、ユニコード文字列を表現できます。SubString型は従来からありましたが、これはString型の一部分を切り出す際に使われています。 将来的にはSubString型もString型に統合される予定のようですが、0.5時点では0.4のASCIIStringUTF8StringString型に統合されたと見るのが良さそうです。

Julia 0.5.0:

julia> typeof("foo")
String

julia> typeof("いえい")
String

Julia 0.4.6:

julia> typeof("foo")
ASCIIString

julia> typeof("いえい")
UTF8String

昔の文字列型はLegacyStrings.jlパッケージに移動され、UTF-16などのエンコーディング方式はStringEncodings.jlパッケージで新たにサポートされるようです。

Fused broadcasting構文

Fused broadcasting構文は、Juliaのベクトル計算を効率化する新しい構文の拡張です。 0.4までのJuliaでは、配列に対して何度も関数適用を行うとその都度新しい配列が生成されていました。 新しいJuliaでは、これが構文レベルで融合(fuse)されるようになります。

sin.(cos.(x))という式があったとしましょう。 0.4ではまずcos.(x)が評価され、xの各要素にcos関数を適用した新しい配列が作られます。 続いてその各要素にsin関数が適用され、また新しい配列が作られ、先ほどの配列は将来的にGCにより破棄されます。 0.5では、この式はまず broadcast(x -> sin(cos(x)), x)のような式に変換されます。 これは、xから値を一つずつ取り出してcossinを順に適用するため、一時的な配列の生成が起きず、一気に2つの関数を適用した配列が生成されます。

さらに、x .= ...のような構文もbroadcast!(identity, x, ...)に変換されるので配列への代入がin-placeに行うことができるようになりました。

マルチスレッドのサポート

マルチスレッドを使った計算の並列化が新たにサポートされました。 今のところ、@threadsマクロを使ってfor文の並列化を行うことができます。 例えば、以下のように複数のタスクを並列に処理することができます。

using Base.Threads

function dotasks(tasks)
    @threads for i in 1:endof(tasks)
        dotask(tasks[i])
    end
end

配列などを領域で分割しておけば配列に対する処理の並列化もできますし、私の作ったBGZFStreams.jlではgzipの解凍処理をこのマルチスレッド機能を使って並列化しています。

テストフレームワークの強化

Julia標準のテストフレームワークがより強化され、テストをきれいに構造化できるようになりました。 今までBase.Testから提供されていた機能では複数のテストをまとめる機能がなく、フラットにすべてのテストを書き連ねるしか方法がありませんでした。 しかし0.5からは@testsetマクロが加わり、次のようにテストを構造化することができるようになりました。

@testset "sum" begin
    @test sum([1,2,3]) == 6
    @test sum([1.0, 2.0, 3.0]) == 6.0
end

その他の変更点

0.5では内部的なものも含め多数の変更がなされています。そのすべてはJuliaのリリースノートから参照することができますので、気になる方は是非調べてみて下さい。

JuliaからRを使う

先日のWACODE夏期講習でRCall.jlのデモをしたら、やはりウケが良かったようなので改めて紹介をします。

RCall.jlはJuliaからR言語の機能を呼び出すツールです。データの受け渡しからREPLでのインタラクティブな実行・プロットも簡単にできます。Juliaを使ってみたいけど、Rの豊富な資産を捨てる訳にはいかないといった方にはピッタリのライブラリです。

インストールは、Juliaの標準的な方法通り、julia -e 'Pkg.update(); Pkg.add(“RCall”)’を実行して下さい。これで最新版のRCall.jlがインストールされることになります。尚、次期Juliaのリリース候補v0.5-RC1では現在動かないようですが、リリースブランチでは直っているのでRC2では使えると思います。

簡単な演算で正しくインストールできたかを確認しましょう。JuliaのREPLを起動して、using RCallとしてRCall.jlを読み込み、R”1 + 1”と打ってR言語1 + 1を実行してみます。以下のように計算結果が表示されればOKです。 f:id:bicycle1885:20160809080928p:plain

R”…"という記法は、Juliaの非標準文字列と呼ばれる機能を利用したものです。この文字列の内側ではRのコードとして解釈され、Rの処理系がコードを実行してくれます。ダブルクウォートが含まれる場合はバックスラッシュでエスケープするか、R””” … “””とトリプルクウォートを使うこともできます。

R側の評価値を取り出す場合には、rcopyが使えます。

julia> rcopy(R"1 + 1")
2.0

julia> rcopy(R"1:5")
5-element Array{Int32,1}:
 1
 2
 3
 4
 5

julia> rcopy(R"""list(x=10, y="foo")""")
Dict{Symbol,Any} with 2 entries:
  :y => "foo"
  :x => 10.0

R側へJuliaの値を渡すには$で変数を埋め込みます。

julia> x = 1
1

julia> R"1 + $x"
RCall.RObject{RCall.RealSxp}
[1] 2

JuliaのDataFrames.jlパッケージが提供しているDataFrameなども自動的にRのdata.frameへと変換してくれます。

julia> using RDatasets

julia> iris = dataset("datasets", "iris");

julia> typeof(iris)
DataFrames.DataFrame

julia> R"head($iris)"
RCall.RObject{RCall.VecSxp}
  SepalLength SepalWidth PetalLength PetalWidth Species
1         5.1        3.5         1.4        0.2  setosa
2         4.9        3.0         1.4        0.2  setosa
3         4.7        3.2         1.3        0.2  setosa
4         4.6        3.1         1.5        0.2  setosa
5         5.0        3.6         1.4        0.2  setosa
6         5.4        3.9         1.7        0.4  setosa

プロットも簡単にできます。以下のようにすると、Juliaから渡したデータをRがプロットしてくれます。

julia> R"""pairs($iris[1:4], pch=21, bg=c("red", "green3", "blue")[unclass($iris$Species)])"""
RCall.RObject{RCall.NilSxp}
NULL

f:id:bicycle1885:20160809081146p:plain

いちいちR”…”で囲んだりが面倒なら、RのREPLに入りましょう。JuliaのREPLから$キーを押すと、RのREPLに入れます。Juliaに戻るにはバックスペースです。 f:id:bicycle1885:20160809081141p:plain

ここではRの文法が自由に使えるので、Rユーザーにとってはやりやすいでしょう。もちろん、Julia側からデータを受け取ることもR”…”の時と同じように可能です。

R> library(ggplot2)

R> ggplot($iris, aes(x=SepalLength, y=SepalWidth, color=Species)) + geom_point()

f:id:bicycle1885:20160809081246p:plain

Rのセッションではヘルプなども参照できますし補完も効きます。Juliaは気になるけどRから離れるのはチョット…という方は是非一度試してみてください。

『R言語徹底解説』(共立出版)をいただきました!

先日、出版社の方から新刊の『R言語徹底解説』(共立出版)をいただきました。ありがとうございます。

R言語徹底解説 / Hadley Wickham 著 石田 基広 市川 太祐 高柳 慎一 福島 真太朗 訳 | 共立出版

f:id:bicycle1885:20160211220632j:plain

この本は、R言語界では知らない人のいないHadley Wickham氏の書いた『Advanced R』の翻訳です。翻訳はこれまた日本のR言語界では知らない人のいないであろう石田先生・市川氏・高柳氏・福島氏の4名の手により行われています。本のページ数は500ページ以上もあり、R言語の仕様・プログラミングスタイル・高速化について書かれています。おそらく、Rの言語自体について書かれた日本語の本としては最も詳細な本でしょう。

プログラミング言語マニアとして見ると、R言語プログラミング言語としては非常に特異な言語だと思います。例えばスカラーにあたる型がなく、すべてベクトルとして扱ったり、引数の遅延評価など、他のメジャーな言語から来た人にとってはとても奇妙な振る舞いをします。また、GNU R (R言語の標準的な処理系) は決して高速とは言えず、実用されるプログラミング言語の中では最も遅い言語の部類でしょう。R言語はどういう動きをするのか、どうしてこうした特徴を持つのかについて答えてくれるのが、この『R言語徹底解説』です。

第6章「関数」では、R言語の関数について詳細に説明しています。6.3節の「すべての操作は関数呼び出しである」では、実はR言語でのオブジェクトに対する操作すべてが関数呼び出しであることを説明しています。ここでは、ifによる分岐やforによる反復などを含む本当にすべての操作が、関数呼び出しにより行われることを強調しています。普通のプログラミング言語では、if文の条件部分の評価のみが先に行われて、その値に応じて2つの別々に処理へと分岐します。関数呼び出しでif文を実現しようと思っても、関数に与えた引数が先に評価されるため、普通は実現できません。しかし、Rでは引数の遅延評価を行うことで、次のように実現可能になるわけです。

> i = 1
> `if`(i == 1, print("yes!"), print("no."))
[1] "yes!"

(6.3節から一部抜粋)

また、未定義の変数を関数に渡すこともできます。

> exists("x")
[1] FALSE
> f <- function(x) {}
> f(x)
NULL

このような機能があるおかげで、例えばggplot2ではggplot(diamonds, aes(x=carat, y=price)) + geom_point()のように簡潔に式が書けるわけです。他にも、この章で説明されている仮引数に対する値の束縛・クロージャ・コピー修正セマンティクスなども、奇妙なRという言語を扱う上で重要な知識でしょう。

第16章では、R言語のパフォーマンスの問題について議論されています。よく知られているようにGNU RによるR言語の実装は速くはありません。これは、R言語がデータ処理のやりやすさに特化した動的な言語であるためです。また、GNU Rの開発は、安定性のためにパフォーマンスを改善する変更に対してはかなり保守的なようです。これは開発において何を優先するかの問題で、決して間違った考えではないと思います。それに、多くの場合についてはR言語のコードを改善するだけでもそれなりにパフォーマンスは良くなります。そのため、第16章から18章のR言語でのパフォーマンスチューニングの手法を学べば十分でしょう。どうしても必要な場合にのみC言語C++を使うことになりますが、これらの言語をR言語から使う方法についても以降の章で説明されています。

個人的な使い方の話をすると、私はパフォーマンスが必要な部分はJulia言語を推しているのですが、R言語の安定感を求める人にはJulia言語はまだ時期尚早かもしれません。私の使い方としては、R言語でデータの基本的な処理 (バイオインフォなのでBioconductorなどのパッケージをよく使う) を行って簡単なテキストファイルを作り、それをもとにJulia言語で計算をして、R言語で可視化するというR > Julia > Rの順で使うことが多いです。なるべくJulia言語でできるところを増やしたいのですが、それでもR言語から離れられるほど成熟していません。

R言語徹底解説』はR言語の基礎を確認してストレスなく使うために役立つ本だと思います。R言語の仕様をあまり意識してこなかった人にとっては、全部は読まなくとも前半部分の第9章くらいまで目を通しておくだけでも随分違うと思います。