Gardner-Altman Plots in Julia

Clement_POIRET · September 29, 2019, 7:18am

Hi everyone,

I’m switching from Python to Julia, and I’m interested in Estimation Statistics to replace classical NHST approach. To do so, I’d like to make Gardner-Altman Plots in Julia but I have to say I’m a bit confused on how to produce it.

Here is a GA Plot:
GAP

So we have two swarmplots representing our 2 independant groups with the effect size and 95% confidence interval on the right plotted over a bootstrap distribution.

In python, we have the library “DABEST”, but I don’t know how to do this in Julia.

Can you help me please?
Thanks,
Clément.

mkborregaard · September 29, 2019, 7:48am

I don’t know of any plotting packages that could currently do this easily out of the box - you’d need to make a function constructing the plot yourself.

nilshg · September 29, 2019, 8:24am

Maybe check out the violin ploys in StatsPlots here for some inspiration - but I agree that you’ll have to do some calculations yourself to construct the bootstrap distribution

baggepinnen · September 29, 2019, 9:21am

Another alternative is to simply call the python code using PyCall

mkborregaard · September 29, 2019, 12:14pm

I would have suggested that as well. Normally you’d just write a recipe if yuo wanted to use (Stats)Plots - but the beeswarm plot PR has not yet been merged on StatsPlots. With beeswarm you could write a userplot recipe with two panels, put the y axis on the right of the second panel, and the legend outside and to the left of the first panel. Then just use group to create the colors. The density plot etc can be found int StatsPlots.

If someone can finish the beeswarm PR I’m happy to add this plot type to StatsPlots

Mattriks · September 29, 2019, 1:35pm

Here’s a Gadfly example, similar to this plot:

using DataFrames, Gadfly, RDatasets, Statistics
iris = dataset("datasets","iris")
irm = by(iris, :Species, :PetalWidth=>mean)
boot(x; n=200) = [mean(rand(x, length(x))) for i in 1:n]
Db = by(iris, :Species, μ=:PetalWidth=>boot)
Db.μdiff = Db.μ - repeat(Db[Db.Species.=="setosa", :μ], outer=3)
labeld = Dict("setosa"=>"", "versicolor"=>"versicolor\n  minus\n setosa", 
  "virginica"=>"virginica \n  minus\n setosa" )
p1 = plot(
  layer(iris, x=:Species, y=:PetalWidth, color=:Species, Geom.beeswarm),
  layer(irm, yintercept=:PetalWidth_mean, Geom.hline(color="grey", style=:dash)),
    Guide.xlabel(nothing), Theme(key_position=:none)
)
p2 = plot(Db, x=:Species, y=:μdiff, color=:Species,
        Geom.violin, Coord.cartesian(ymin=1, ymax=2),
    Scale.x_discrete(labels=i->labeld[i]),
    Guide.xlabel(nothing), Guide.ylabel("Mean difference"), 
    Theme(key_position=:none)
)
vstack(p1, p2)