using RCall, DataFrames

@rimport base as R
R.Sys_setenv(LANG = "en") # こうしておかないとR側でエラーが出たときにハングする

# 以下はJulia側でggplot2を便利に使うための設定

rcalljl_options(; kwargs...) = rcopy(RCall.rcall_p(:options, rcalljl_options=Dict(kwargs)))

function rplotsize(w, h)
    if RCall.ijulia_mime == MIME("image/svg+xml")
        rcalljl_options(; width=w/100, height=h/100)
    else
        rcalljl_options(; width=w, height=h)
    end
end

function rplotpng(; kwargs...)
    RCall.ijulia_setdevice(MIME("image/png"); kwargs...)
    RCall.ijulia_mime, rcalljl_options()
end

function rplotsvg(; kwargs...)
    RCall.ijulia_setdevice(MIME("image/svg+xml"); kwargs...)
    RCall.ijulia_mime, rcalljl_options()
end

# Avoid svg display bug
using Base64
function RCall.ijulia_displayfile(m::MIME"image/svg+xml", f)
    open(f) do io
        svg = read(io, String)
        base64 = base64encode(svg)
        html = """<img src="data:$m;base64,$base64" />"""
        display(MIME("text/html"), html)
    end
end

@rlibrary ggplot2
rplotsvg()

rplotsize(640, 400) # グラフの表示サイズを設定
t = range(-10, 10, length=1000)
data = DataFrame(t = t, x = cos.(t), y = sin.(t))
ggplot(data=data, aes(x=:t, y=:x)) + 
geom_line(color="red") +
geom_line(aes(y=:y), color="blue")

@rlibrary stats # これでRのfisher_testなどをJuliaで使用可能

A = [
    10 10
     7 27
]

fisher_result_R = fisher_test(A) # P値が5%未満なのに, 95%信頼区間が1を含む

fisher_result = rcopy(fisher_result_R) # RのデータをJuliaのデータに変換9

# fisher_testの自前での実装 (エラーチェックの類に欠けた手抜きであることに注意！)

module My # 使い捨てモジュール

using Distributions, Roots
⪅(x ,y) = x < y || x ≈ y

or(a, b, c, d) = a*d/(b*c)
or(x::AbstractVecOrMat) = or(x...)

pval(d::FisherNoncentralHypergeometric, k) =
    sum(pdf(d, j) for j in support(d) if pdf(d,j) ⪅ pdf(d, k))
pval(a, b, c, d, ω=1.0) =
    pval(FisherNoncentralHypergeometric(a+b, c+d, a+c, ω), a)
pval(A::AbstractVecOrMat, ω=1.0) = pval(A..., ω)

ci(a, b, c, d, α=0.05) = 
    exp.(find_zeros(t -> pval(a, b, c, d, exp(t)) - α, -1e2, 1e2))
ci(A::AbstractVecOrMat, α=0.05) = ci(A..., α)

struct FisherTest{D, T}
    data::D
    ω::T
    α::T
    odds_ratio::T
    p_value::T
    conf_int::Vector{T}
end

Base.show(io::IO, x::FisherTest) = print(io, """
            Fisher's Exact Test for Count Data
    
    data: $(x.data)
    p-value: $(x.p_value)
    alternative hypothesis: true odds ratio is not equal to $(x.ω)
    $(100(1 - x.α)) percent confidence interval:
        $(x.conf_int)
    odds ratio: $(x.odds_ratio)
    """)

function fisher_test(data; ω=1.0, α = 0.05)
    odds_ratio = or(data)
    p_value = pval(data, ω)
    conf_int = ci(data, α)
    FisherTest(data, ω, α, odds_ratio, p_value, conf_int)
end

end

My_fisher_result = My.fisher_test(A)

dump(My_fisher_result)