Skip to content

Commit

Permalink
more plot recipes
Browse files Browse the repository at this point in the history
  • Loading branch information
alexpkeil1 committed Sep 15, 2023
1 parent aceab7e commit 5fadd17
Show file tree
Hide file tree
Showing 2 changed files with 97 additions and 5 deletions.
14 changes: 12 additions & 2 deletions src/examples.jl
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,6 @@

using LSurvival, LinearAlgebra, RCall, BenchmarkTools, Random
using Plots
using RecipesBase

###################################################################
# Test data
Expand Down Expand Up @@ -39,9 +38,20 @@ dat4 = ( # mine
###################################################################
# plotting survival data
###################################################################
id, int, outt, data =
LSurvival.dgm(MersenneTwister(123123), 100, 100; afun = LSurvival.int_0)
data[:, 1] = round.(data[:, 1], digits = 3)
d, X = data[:, 4], data[:, 1:3]
wt = rand(length(d))
wt ./= (sum(wt) / length(wt))

# equivalent methods for unweighted, default efron partial likelihood
ft = coxph(X, int, outt, d, id=ID.(id))
k = kaplan_meier(int, outt, d, id=ID.(id))


plot(R, maxids=10)
plot(ft.R, maxids=10)
plot(k)

###################################################################
# Fitting a basic Cox model, Kaplan-Meier curve using preferred functions
Expand Down
88 changes: 85 additions & 3 deletions src/plot_recipes.jl
Original file line number Diff line number Diff line change
Expand Up @@ -66,7 +66,7 @@ R = LSurvivalResp(dat4.enter, dat4.exit, dat4.status)
```
"""
@recipe function f(r::KMSurv)
@recipe function f(k::KMSurv)
# global
xlabel --> "Time" # --> sets default
ylabel --> "Survival"
Expand Down Expand Up @@ -163,7 +163,7 @@ end


# baseline hazard plot for Cox model

@userplot BaseHazPlot
"""
Plotting baseline hazard for a Cox model
Expand All @@ -183,7 +183,8 @@ plot!(ftb, label="Breslow")
```
"""
@recipe function f(ft::PHModel)
@recipe function f(h::BaseHazPlot)
ft = h.args[1]
# global
xlabel --> "Time" # --> sets default
ylabel --> "Cumulative hazard"
Expand All @@ -202,3 +203,84 @@ plot!(ftb, label="Breslow")
vcat(minT, times, maxT), vcat(0, cumbasehaz, maximum(cumbasehaz))
end
end



@userplot CoxDX
"""
```julia
using Plots, LSurvival
dat2 = (
enter = [1, 2, 5, 2, 1, 7, 3, 4, 8, 8],
exit = [2, 3, 6, 7, 8, 9, 9, 9, 14, 17],
status = [1, 1, 1, 1, 1, 1, 1, 0, 0, 0],
x = [1, 0, 0, 1, 0, 1, 1, 1, 0, 0],
)
fte = coxph(@formula(Surv(enter, exit, status)~x), dat2)
coxdx(fte)
```
"""
@recipe function f(h::CoxDX; par=1)
ft = h.args[1]
time = ft.bh[:,1]
nms = coefnames(ft)
xlab --> "Time"
ylab --> "Schoenfeld residual"
res = residuals(ft, type = "schoenfeld")
@series begin
seriestype := :scatter
shape := :circle
markerstrokewidth := 0
color --> :black
label --> nms[par]
markesize --> 2
markeralpha --> 0.5
time, res[:,par]
end
end

@userplot CoxInfluence
"""
```julia
using Plots, LSurvival
dat2 = (
enter = [1, 2, 5, 2, 1, 7, 3, 4, 8, 8],
exit = [2, 3, 6, 7, 8, 9, 9, 9, 14, 17],
status = [1, 1, 1, 1, 1, 1, 1, 0, 0, 0],
x = [1, 0, 0, 1, 0, 1, 1, 1, 0, 0],
)
fte = coxph(@formula(Surv(enter, exit, status)~x), dat2)
coxinfluence(fte, type="jackknife")
```
"""
@recipe function f(h::CoxInfluence; type="dfbeta")
!issubset([type], ["dfbeta", "dfbetas", "jackknife"]) && throw("type must be 'dfbeta', 'dfbetas' or 'jackknife'")
ft = h.args[1]
id = values(ft.R.id)
nms = coefnames(ft)
xlab --> "ID"
ylab --> "Residual"
res = residuals(ft, type = type)
grid --> false
@series begin
seriestype := :scatter
shape := :circle
markerstrokewidth := 0
color --> :black
label --> type
markesize --> 2
markeralpha --> 0.5
id, res
end
@series begin
seriestype := :hline
label := ""
color := :black
[0]
end
end

0 comments on commit 5fadd17

Please sign in to comment.