Kapitel 80 Ordinale Logistische Regression
Dieses Kapitel ist noch nicht fertig und ist unkorrigiert.
80.0.1 Wir benötigen folgende Pakete:
library(dplyr)
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
library(ggplot2)
library(Statamarkdown)
## Stata found at C:/Program Files/Stata17/StataMP-64.exe
## The 'stata' engine is ready to use.
library(epiR)
## Loading required package: survival
## Package epiR 2.0.53 is loaded
## Type help(epi.about) for summary information
## Type browseVignettes(package = 'epiR') to learn how to use epiR for applied epidemiological analyses
##
library(sjPlot)
library(rms)
## Loading required package: Hmisc
## Loading required package: lattice
## Loading required package: Formula
##
## Attaching package: 'Hmisc'
## The following objects are masked from 'package:dplyr':
##
## src, summarize
## The following objects are masked from 'package:base':
##
## format.pval, units
## Loading required package: SparseM
##
## Attaching package: 'SparseM'
## The following object is masked from 'package:base':
##
## backsolve
library(summarytools)
##
## Attaching package: 'summarytools'
## The following objects are masked from 'package:Hmisc':
##
## label, label<-
80.1 Daten simulieren
Für unser erstes Beispiel simulieren wir Daten. Um das Verständnis zu fördern, simulieren wir zuerst zwei binäre Variablen.
=666
set.seed=1000
n <-1:100
id<-sample(c(0,1), size=n, rep=TRUE, prob=c(0.6, 0.4))
osteoarthritis
<-data.frame(id, osteoarthritis)
data
<-data %>%
datamutate(pain_yes_no=case_when(
==1~sample(c(0,1), size=length(osteoarthritis), rep=TRUE, prob=c(0.8, 0.2)),
osteoarthritis==0~sample(c(0,1), size=length(osteoarthritis), rep=TRUE, prob=c(0.95, 0.05))))
osteoarthritis::export(data, "bin_log_ordinal.dta") rio
Schauen wir uns zuerst eine logistische Regression an:
<-glm(pain_yes_no~osteoarthritis, data=data, family='binomial')
logregtab_model(logreg)
pain yes no | |||
---|---|---|---|
Predictors | Odds Ratios | CI | p |
(Intercept) | 0.06 | 0.04 – 0.08 | <0.001 |
osteoarthritis | 4.79 | 3.16 – 7.44 | <0.001 |
Observations | 1000 | ||
R2 Tjur | 0.059 |
Nun eine ordinale logistische Regression (die natürlich hier keinen Sinn macht, aber das hilft uns beim Verstehen). Der Befehl polr aus dem MASS Paket funktioniert nur, wenn die abhängige Variable mehr als zwei Kategorien hat.
<-MASS::polr(factor(pain_yes_no)~osteoarthritis, data=data) or.fit
Error in MASS::polr(factor(pain_yes_no) ~ osteoarthritis, data = data): response must have 3 or more levels
tab_model(or.fit)
Error in tab_model(or.fit): object 'or.fit' not found
Wir können den orm Befehl aus dem rms Paket von Frank Harrell benutzen.
<-rms::orm(pain_yes_no~osteoarthritis, data=data)
orm.fittab_model(orm.fit)
pain yes no | |||
---|---|---|---|
Predictors | Odds Ratios | CI | p |
Intercept | 0.06 | 0.04 – 0.08 | <0.001 |
osteoarthritis | 4.79 | 3.12 – 7.34 | <0.001 |
Observations | 1000 | ||
R2 | 1.000 |
Oder den Befehl ologit in Stata:
use bin_log_ordinal.dta
ologit pain_yes_no osteoarthritis ,or
Iteration 0: log likelihood = -366.92499
Iteration 1: log likelihood = -339.8507
Iteration 2: log likelihood = -337.41003
Iteration 3: log likelihood = -337.40342
Iteration 4: log likelihood = -337.40342
Ordered logistic regression Number of obs = 1,000
LR chi2(1) = 59.04
Prob > chi2 = 0.0000
Log likelihood = -337.40342 Pseudo R2 = 0.0805
------------------------------------------------------------------------------
pain_yes_no | Odds ratio Std. err. z P>|z| [95% conf. interval]
-------------+----------------------------------------------------------------
osteoarthr~s | 4.788941 1.043921 7.19 0.000 3.123846 7.341576
-------------+----------------------------------------------------------------
/cut1 | 2.860414 .1817661 2.504159 3.216669
------------------------------------------------------------------------------
Note: Estimates are transformed only in the first equation to odds ratios.
Wir sehen, dass die Odds Ratio bei beiden Modellen praktisch gleich ist.
Wir simulieren etwas komplexere Daten.
=666
set.seed=1000
n<-1:n
id<-sample(c(0,1,2), size=n, rep=TRUE, prob=c(0.6, 0.4,0.2))
Anzahl_Gelenke_Mit_Arthrose
<-data.frame(id, Anzahl_Gelenke_Mit_Arthrose)
data
<-data %>%
datamutate(pain_no_sometimes_often=case_when(
==2~sample(c(0,1,2), size=length(Anzahl_Gelenke_Mit_Arthrose), rep=TRUE, prob=c(0.6, 0.25, 0.15)),
Anzahl_Gelenke_Mit_Arthrose==1~sample(c(0,1,2), size=length(Anzahl_Gelenke_Mit_Arthrose), rep=TRUE, prob=c(0.7, 0.2,0.1)),
Anzahl_Gelenke_Mit_Arthrose==0~sample(c(0,1,2), size=length(Anzahl_Gelenke_Mit_Arthrose), rep=TRUE, prob=c(0.9, 0.05, 0.025))))
Anzahl_Gelenke_Mit_Arthrose::export(data, "ordinal_log_ordinal.dta") rio
Zuerst wieder in R mit dem polr Befehl:
<-MASS::polr(factor(pain_no_sometimes_often)~factor(Anzahl_Gelenke_Mit_Arthrose), data=data)
polr.fittab_model(polr.fit)
Re-fitting to get Hessian
factor(pain no sometimes often) |
|||
---|---|---|---|
Predictors | Odds Ratios | CI | p |
0|1 | 14.79 | 10.28 – 21.29 | <0.001 |
1|2 | 55.96 | 36.78 – 85.14 | <0.001 |
Anzahl Gelenke Mit Arthrose [1] |
6.54 | 4.30 – 10.21 | <0.001 |
Anzahl Gelenke Mit Arthrose [2] |
11.99 | 7.63 – 19.27 | <0.001 |
Observations | 1000 | ||
R2 Nagelkerke | 0.192 |
Und mit dem orm Befehl.
library(rms)
<-data %>%
datamutate(Anzahl_Gelenke_Mit_Arthrose=factor(Anzahl_Gelenke_Mit_Arthrose))
<-rms::orm(pain_no_sometimes_often~Anzahl_Gelenke_Mit_Arthrose, data=data)
orm.fittab_model(orm.fit)
pain no sometimes often | |||
---|---|---|---|
Predictors | Odds Ratios | CI | p |
y>=1 | 0.07 | 0.05 – 0.10 | <0.001 |
Anzahl Gelenke Mit Arthrose=1 |
6.54 | 4.25 – 10.07 | <0.001 |
Anzahl Gelenke Mit Arthrose=2 |
11.99 | 7.55 – 19.03 | <0.001 |
Observations | 1000 | ||
R2 | 1.000 |
Und wieder mit Stata:
use ordinal_log_ordinal.dta
ologit pain_no_sometimes_often i.Anzahl_Gelenke_Mit_Arthrose ,or
Iteration 0: log likelihood = -656.08703
Iteration 1: log likelihood = -587.72574
Iteration 2: log likelihood = -580.52438
Iteration 3: log likelihood = -580.4747
Iteration 4: log likelihood = -580.4747
Ordered logistic regression Number of obs = 1,000
LR chi2(2) = 151.22
Prob > chi2 = 0.0000
Log likelihood = -580.4747 Pseudo R2 = 0.1152
------------------------------------------------------------------------------
pain_no_so~n | Odds ratio Std. err. z P>|z| [95% conf. interval]
-------------+----------------------------------------------------------------
Anzahl_Gel~e |
1 | 6.539586 1.438919 8.53 0.000 4.248725 10.06565
2 | 11.99118 2.826658 10.54 0.000 7.554571 19.0333
-------------+----------------------------------------------------------------
/cut1 | 2.694182 .1855077 2.330593 3.05777
/cut2 | 4.024605 .2138546 3.605457 4.443752
------------------------------------------------------------------------------
Note: Estimates are transformed only in the first equation to odds ratios.
Mit Stata können wir den Koeffizienten auch als logOdds ausgeben lassen:
use ordinal_log_ordinal.dta
ologit pain_no_sometimes_often i.Anzahl_Gelenke_Mit_Arthrose
Iteration 0: log likelihood = -656.08703
Iteration 1: log likelihood = -587.72574
Iteration 2: log likelihood = -580.52438
Iteration 3: log likelihood = -580.4747
Iteration 4: log likelihood = -580.4747
Ordered logistic regression Number of obs = 1,000
LR chi2(2) = 151.22
Prob > chi2 = 0.0000
Log likelihood = -580.4747 Pseudo R2 = 0.1152
------------------------------------------------------------------------------
pain_no_so~n | Coefficient Std. err. z P>|z| [95% conf. interval]
-------------+----------------------------------------------------------------
Anzahl_Gel~e |
1 | 1.877874 .2200321 8.53 0.000 1.446619 2.309129
2 | 2.484171 .2357281 10.54 0.000 2.022153 2.94619
-------------+----------------------------------------------------------------
/cut1 | 2.694182 .1855077 2.330593 3.05777
/cut2 | 4.024605 .2138546 3.605457 4.443752
------------------------------------------------------------------------------
In Stata können wir auch die Wahrscheinlichkeiten für jede Kategorie schätzen, hier für die Antwortskategorie “oft” Schmerzen (pain_no_sometimes_often == 2).
use ordinal_log_ordinal.dta
ologit pain_no_sometimes_often i.Anzahl_Gelenke_Mit_Arthrose
at(Anzahl_Gelenke_Mit_Arthrose=(0(1)2)) predict(outcome(2)) atmeans margins,
Iteration 0: log likelihood = -656.08703
Iteration 1: log likelihood = -587.72574
Iteration 2: log likelihood = -580.52438
Iteration 3: log likelihood = -580.4747
Iteration 4: log likelihood = -580.4747
Ordered logistic regression Number of obs = 1,000
LR chi2(2) = 151.22
Prob > chi2 = 0.0000
Log likelihood = -580.4747 Pseudo R2 = 0.1152
------------------------------------------------------------------------------
pain_no_so~n | Coefficient Std. err. z P>|z| [95% conf. interval]
-------------+----------------------------------------------------------------
Anzahl_Gel~e |
1 | 1.877874 .2200321 8.53 0.000 1.446619 2.309129
2 | 2.484171 .2357281 10.54 0.000 2.022153 2.94619
-------------+----------------------------------------------------------------
/cut1 | 2.694182 .1855077 2.330593 3.05777
/cut2 | 4.024605 .2138546 3.605457 4.443752
------------------------------------------------------------------------------
Adjusted predictions Number of obs = 1,000
Model VCE: OIM
Expression: Pr(pain_no_sometimes_often==2), predict(outcome(2))
1._at: Anzahl_Gelenke_Mit_Arthrose = 0
2._at: Anzahl_Gelenke_Mit_Arthrose = 1
3._at: Anzahl_Gelenke_Mit_Arthrose = 2
------------------------------------------------------------------------------
| Delta-method
| Margin std. err. z P>|z| [95% conf. interval]
-------------+----------------------------------------------------------------
_at |
1 | .0175567 .0036887 4.76 0.000 .0103271 .0247864
2 | .1046371 .0141246 7.41 0.000 .0769535 .1323207
3 | .1764723 .0240215 7.35 0.000 .1293911 .2235535
------------------------------------------------------------------------------
wir können die geschätzten Wahrscheinlichkeiten noch mit den beobachteten Häufigkeiten verlgeichen.
use ordinal_log_ordinal.dta
tab Anzahl_Gelenke_Mit_Arthrose pain_no_sometimes_often, row
| Key |
|----------------|
| frequency |
| row percentage |
+----------------+
Anzahl_Gel |
enke_Mit_A | pain_no_sometimes_often
rthrose | 0 1 2 | Total
-----------+---------------------------------+----------
0 | 458 23 8 | 489
| 93.66 4.70 1.64 | 100.00
-----------+---------------------------------+----------
1 | 228 67 34 | 329
| 69.30 20.36 10.33 | 100.00
-----------+---------------------------------+----------
2 | 101 48 33 | 182
| 55.49 26.37 18.13 | 100.00
-----------+---------------------------------+----------
Total | 787 138 75 | 1,000
| 78.70 13.80 7.50 | 100.00
80.2 Weitere Beispiele
Wir simulieren weitere Daten.
set.seed
[1] 666
<-sample(c(0,1), size=1000, rep=TRUE, prob=c(0.6,0.4))
osteoarthritis<-1:length(osteoarthritis)
id<-data.frame(id, osteoarthritis)
data
<-data %>%
datamutate(walking_difficulty=factor(case_when(
==0~sample(c("not difficult","quite difficult","very difficult","impossible"), size=length(osteoarthritis), rep=TRUE, prob=c(0.8, 0.1, 0.025, 0.025)),
osteoarthritis==1~sample(c("not difficult","quite difficult","very difficult","impossible"), size=length(osteoarthritis), rep=TRUE, prob=c(0.725, 0.15, 0.075, 0.05))))) %>%
osteoarthritismutate(walking_difficulty=factor(walking_difficulty, levels=c("not difficult","quite difficult","very difficult","impossible")))
levels(data$walking_difficulty)
[1] "not difficult" "quite difficult" "very difficult" "impossible"
Hier die Verteilung der Daten:
ctable(factor(data$osteoarthritis), data$walking_difficulty)
Cross-Tabulation, Row Proportions
factor(data$osteoarthritis) * walking_difficulty
----------------------------- -------------------- --------------- ----------------- ---------------- ------------ ---------------
walking_difficulty not difficult quite difficult very difficult impossible Total
factor(data$osteoarthritis)
0 517 (85.3%) 63 (10.4%) 12 (2.0%) 14 (2.3%) 606 (100.0%)
1 295 (74.9%) 53 (13.5%) 23 (5.8%) 23 (5.8%) 394 (100.0%)
Total 812 (81.2%) 116 (11.6%) 35 (3.5%) 37 (3.7%) 1000 (100.0%)
----------------------------- -------------------- --------------- ----------------- ---------------- ------------ ---------------
::export(data, "oa_walking.dta") rio
use oa_walking.dta
ologit walking_difficulty i.osteoarthritis
Iteration 0: log likelihood = -658.3034
Iteration 1: log likelihood = -649.12874
Iteration 2: log likelihood = -648.99698
Iteration 3: log likelihood = -648.99696
Ordered logistic regression Number of obs = 1,000
LR chi2(1) = 18.61
Prob > chi2 = 0.0000
Log likelihood = -648.99696 Pseudo R2 = 0.0141
------------------------------------------------------------------------------
walking_di~y | Coefficient Std. err. z P>|z| [95% conf. interval]
-------------+----------------------------------------------------------------
1.osteoart~s | .7009046 .1624187 4.32 0.000 .3825697 1.019239
-------------+----------------------------------------------------------------
/cut1 | 1.771396 .1144781 1.547023 1.995768
/cut2 | 2.879737 .1489003 2.587898 3.171576
/cut3 | 3.589341 .188595 3.219702 3.958981
------------------------------------------------------------------------------
use oa_walking.dta
ologit walking_difficulty i.osteoarthritis, or
Iteration 0: log likelihood = -658.3034
Iteration 1: log likelihood = -649.12874
Iteration 2: log likelihood = -648.99698
Iteration 3: log likelihood = -648.99696
Ordered logistic regression Number of obs = 1,000
LR chi2(1) = 18.61
Prob > chi2 = 0.0000
Log likelihood = -648.99696 Pseudo R2 = 0.0141
------------------------------------------------------------------------------
walking_di~y | Odds ratio Std. err. z P>|z| [95% conf. interval]
-------------+----------------------------------------------------------------
1.osteoart~s | 2.015575 .3273671 4.32 0.000 1.466047 2.771086
-------------+----------------------------------------------------------------
/cut1 | 1.771396 .1144781 1.547023 1.995768
/cut2 | 2.879737 .1489003 2.587898 3.171576
/cut3 | 3.589341 .188595 3.219702 3.958981
------------------------------------------------------------------------------
Note: Estimates are transformed only in the first equation to odds ratios.
use oa_walking.dta
tab walking_difficulty
tab osteoarthritis walking_difficulty, nol col
ologit walking_difficulty i.osteoarthritis, or
walking_difficu |
lty | Freq. Percent Cum.
----------------+-----------------------------------
not difficult | 812 81.20 81.20
quite difficult | 116 11.60 92.80
very difficult | 35 3.50 96.30
impossible | 37 3.70 100.00
----------------+-----------------------------------
Total | 1,000 100.00
+-------------------+
| Key |
|-------------------|
| frequency |
| column percentage |
+-------------------+
osteoarthr | walking_difficulty
itis | 1 2 3 4 | Total
-----------+--------------------------------------------+----------
0 | 517 63 12 14 | 606
| 63.67 54.31 34.29 37.84 | 60.60
-----------+--------------------------------------------+----------
1 | 295 53 23 23 | 394
| 36.33 45.69 65.71 62.16 | 39.40
-----------+--------------------------------------------+----------
Total | 812 116 35 37 | 1,000
| 100.00 100.00 100.00 100.00 | 100.00
Iteration 0: log likelihood = -658.3034
Iteration 1: log likelihood = -649.12874
Iteration 2: log likelihood = -648.99698
Iteration 3: log likelihood = -648.99696
Ordered logistic regression Number of obs = 1,000
LR chi2(1) = 18.61
Prob > chi2 = 0.0000
Log likelihood = -648.99696 Pseudo R2 = 0.0141
------------------------------------------------------------------------------
walking_di~y | Odds ratio Std. err. z P>|z| [95% conf. interval]
-------------+----------------------------------------------------------------
1.osteoart~s | 2.015575 .3273671 4.32 0.000 1.466047 2.771086
-------------+----------------------------------------------------------------
/cut1 | 1.771396 .1144781 1.547023 1.995768
/cut2 | 2.879737 .1489003 2.587898 3.171576
/cut3 | 3.589341 .188595 3.219702 3.958981
------------------------------------------------------------------------------
Note: Estimates are transformed only in the first equation to odds ratios.
80.2.1 Diagnostic of the proportional odds assumption
Der erste Test können wir mit dem Befehl omodel durchführen. omodel müssen wir zuerst mit ssc install omodel instalieren. Wir sehen, dass der Test nicht statististisch signfikant ist, das bedeutet, dass die Voraussetzungen für die ordinale Regression gegeben ist.
use oa_walking.dta
logit walking_difficulty osteoarthritis omodel
Iteration 0: log likelihood = -658.3034
Iteration 1: log likelihood = -649.12874
Iteration 2: log likelihood = -648.99698
Iteration 3: log likelihood = -648.99696
Ordered logit estimates Number of obs = 1000
LR chi2(1) = 18.61
Prob > chi2 = 0.0000
Log likelihood = -648.99696 Pseudo R2 = 0.0141
------------------------------------------------------------------------------
walking_di~y | Coefficient Std. err. z P>|z| [95% conf. interval]
-------------+----------------------------------------------------------------
osteoarthr~s | .7009045 .1624187 4.32 0.000 .3825696 1.019239
-------------+----------------------------------------------------------------
_cut1 | 1.771396 .1144781 (Ancillary parameters)
_cut2 | 2.879737 .1489003
_cut3 | 3.589342 .188595
------------------------------------------------------------------------------
Approximate likelihood-ratio test of proportionality of odds
across response categories:
chi2(2) = 4.66
Prob > chi2 = 0.0973
Für den zweiten Test müssen wir zuerst spost installieren (search spost) Auch dieser Test ist nihct signifikant und gibt uns an, dass wir mit ologit die Voraussetzungen nicht verletzen, dass die Odds Ratio gleich ist für jede Antwortskategorie.
use oa_walking.dta
ologit walking_difficulty i.osteoarthritis
brant
Iteration 0: log likelihood = -658.3034
Iteration 1: log likelihood = -649.12874
Iteration 2: log likelihood = -648.99698
Iteration 3: log likelihood = -648.99696
Ordered logistic regression Number of obs = 1,000
LR chi2(1) = 18.61
Prob > chi2 = 0.0000
Log likelihood = -648.99696 Pseudo R2 = 0.0141
------------------------------------------------------------------------------
walking_di~y | Coefficient Std. err. z P>|z| [95% conf. interval]
-------------+----------------------------------------------------------------
1.osteoart~s | .7009046 .1624187 4.32 0.000 .3825697 1.019239
-------------+----------------------------------------------------------------
/cut1 | 1.771396 .1144781 1.547023 1.995768
/cut2 | 2.879737 .1489003 2.587898 3.171576
/cut3 | 3.589341 .188595 3.219702 3.958981
------------------------------------------------------------------------------
Brant test of parallel regression assumption
| chi2 p>chi2 df
------------------+------------------------------
All | 3.98 0.137 2
------------------+------------------------------
1.osteoarthritis | 3.98 0.137 2
A significant test statistic provides evidence that the parallel
regression assumption has been violated.
Wir müssten hier eigentlich eine nomimale logistische Regression durchführen.
use oa_walking.dta
mlogit walking_difficulty i.osteoarthritis
Iteration 0: log likelihood = -658.3034
Iteration 1: log likelihood = -647.47231
Iteration 2: log likelihood = -646.91571
Iteration 3: log likelihood = -646.91426
Iteration 4: log likelihood = -646.91426
Multinomial logistic regression Number of obs = 1,000
LR chi2(3) = 22.78
Prob > chi2 = 0.0000
Log likelihood = -646.91426 Pseudo R2 = 0.0173
------------------------------------------------------------------------------
walking_di~y | Coefficient Std. err. z P>|z| [95% conf. interval]
-------------+----------------------------------------------------------------
not_diffic~t | (base outcome)
-------------+----------------------------------------------------------------
quite_diff~t |
1.osteoart~s | .3882247 .2001625 1.94 0.052 -.0040865 .7805359
_cons | -2.104908 .1334438 -15.77 0.000 -2.366453 -1.843363
-------------+----------------------------------------------------------------
very_diffi~t |
1.osteoart~s | 1.211655 .3635047 3.33 0.001 .499199 1.924111
_cons | -3.763136 .2920061 -12.89 0.000 -4.335458 -3.190815
-------------+----------------------------------------------------------------
impossible |
1.osteoart~s | 1.057504 .3467433 3.05 0.002 .3779001 1.737109
_cons | -3.608986 .2708557 -13.32 0.000 -4.139853 -3.078118
------------------------------------------------------------------------------
Wir simulieren neue Daten:
=666
set.seed<-sample(c(0,1), size=1000, rep=TRUE, prob=c(0.6,0.4))
osteoarthritis<-1:length(osteoarthritis)
id<-data.frame(id, osteoarthritis)
data
<-data %>%
datamutate(walking_difficulty=factor(case_when(
==0~sample(c("not difficult","quite difficult","very difficult","impossible"), size=length(osteoarthritis), rep=TRUE, prob=c(0.6, 0.1, 0.2, 0.1)),
osteoarthritis==1~sample(c("not difficult","quite difficult","very difficult","impossible"), size=length(osteoarthritis), rep=TRUE, prob=c(0.40, 0.15, 0.3, 0.15))))) %>%
osteoarthritismutate(walking_difficulty=factor(walking_difficulty, levels=c("not difficult","quite difficult","very difficult","impossible")))
::export(data, "oa_walking.dta") rio
use oa_walking.dta
ologit walking_difficulty i.osteoarthritis
brant
Iteration 0: log likelihood = -1237.074
Iteration 1: log likelihood = -1208.3593
Iteration 2: log likelihood = -1208.2851
Iteration 3: log likelihood = -1208.2851
Ordered logistic regression Number of obs = 1,000
LR chi2(1) = 57.58
Prob > chi2 = 0.0000
Log likelihood = -1208.2851 Pseudo R2 = 0.0233
------------------------------------------------------------------------------
walking_di~y | Coefficient Std. err. z P>|z| [95% conf. interval]
-------------+----------------------------------------------------------------
1.osteoart~s | .9139844 .121393 7.53 0.000 .6760586 1.15191
-------------+----------------------------------------------------------------
/cut1 | .2519758 .0809125 .0933902 .4105614
/cut2 | .7485291 .0840961 .5837037 .9133545
/cut3 | 2.215203 .1104465 1.998732 2.431674
------------------------------------------------------------------------------
Brant test of parallel regression assumption
| chi2 p>chi2 df
------------------+------------------------------
All | 0.07 0.964 2
------------------+------------------------------
1.osteoarthritis | 0.07 0.964 2
A significant test statistic provides evidence that the parallel
regression assumption has been violated.
<-rio::import("https://doi.org/10.1371/journal.pone.0270030.s001") data
Error: Format not supported