Menurut Agresti (2002) model regresi logistik termasuk dalam model linear umum (Generalized Linear Models). Model regresi logistik juga dapat disebut sebagai model logit. Model logit digunakan untuk memodelkan hubungan antara satu variabel respon yang bersifat kategori dan beberapa variabel bebas yang bersifat kategori maupun kontinu. Apabila variabel respon terdiri dari lebih dari dua kategori dan terdapat tingkatan dalam kategori tersebut (skala ordinal) maka dinamakan model regresi logistik ordinal. Dalam Agresti (2002) model untuk regresi logistik ordinal adalah model logit kumulatif (cumulative logit models). Pada model logit ini sifat ordinal dari respon Y dituangkan dalam peluang kumulatif. Terdapat dua tujuan utama dalam regresi logistik, yaitu mengidentifikasi variabel independen yang mempengaruhi pengelompokkan variabel dependen dan membuat sistem klasifikasi yang berbasis pada model logistik untuk pengelompokkan.
Misalkan variabel respon Y memiliki J buah kategori berskala ordinal dan \(x_i\) menyatakan vektor variabel prediktor ke-p pada pengamatan ke-i, \(\pmb{x_i}=\begin{bmatrix} x_{i1} & x_{i2} & ... & x_{ip} \end{bmatrix}^T\) dengan \(i=1,2,...,n\), maka model logit kumulatif dinyatakan sebagai berikut.
\[ \boxed{logit[P(Y_i\leq j|\pmb{x_i})]=\alpha_j+\pmb{x_i}^T\pmb{\beta}} \]
dengan \(P(Y_i\leq j|\pmb{x_i})\) adalah peluang kumulatif kurang dari atau sama dengan kategori ke-j apabila diketahui \(x_i\), \(\alpha_j\) merupakan parameter intersep yang meningkat seiring j \(\alpha_1\leq \alpha_2 \leq ... \leq \alpha_{J-1}\) karena \(P(Y_i\leq j|\pmb{x_i})\) meningkat seiring j untuk \(x\) yang tetap, dan \(\pmb{\beta}=\begin{bmatrix} \beta_1 & \beta_2 & ... &\beta_p \end{bmatrix}^T\) merupakan vektor koefisien regresi yang bersesuaian dengan \(x_1, x_2, ..., x_p\).
Logit kumulatif didefinisikan sebagai \[ logit[P(Y_i\leq j|\pmb{x_i})]=ln \begin{bmatrix} \frac{P(Y_i\leq j|\pmb{x_i})} {1-P(Y_i\leq j|\pmb{x_i})} \end{bmatrix} , \quad g=1,2,...,J-1 \]
maka model regresi logistik ordinal dapat dinyatakan sebagai \[ logit[P(Y_i\leq j|\pmb{x_i})]=ln \begin{bmatrix} \frac{P(Y_i\leq j|\pmb{x_i})} {1-P(Y_i\leq j|\pmb{x_i})} \end{bmatrix} = \alpha_j+\pmb{x_i}^T\pmb{\beta} \]
sehingga diperoleh \[ P(Y_i\leq j|\pmb{x_i})= \frac{exp(\alpha_j+\pmb{x_i}^T\pmb{\beta})} {1+exp(\alpha_j+\pmb{x_i}^T\pmb{\beta})} , \quad j=1,2,...,J-1 \]
Misalkan \(\pi_j(\pmb{x_i})=P(Y_i\leq j|\pmb{x})\) menyatakan peluang variabel respon pada pengamatan ke-i kategori ke-j jika diketahui \(\pmb{x_i}\), maka \[ \begin{aligned} P(Y_i\leq j|\pmb{x_i})&= P(Y_i=1|\pmb{x_i})+P(Y_i=2|\pmb{x_i})+...+P(Y_i=j|\pmb{x_i}) \\ &=\pi_1(\pmb{x_i})+\pi_2(\pmb{x_i})+...+\pi_j(\pmb{x_i}) \end{aligned} \]
sehingga peluang untuk masing-masing ketagori respon dapat dinyatakan sebagai \[ \pi_j(\pmb{x_i})=P(Y_i=j|\pmb{x_i})-P(Y_i=j-1|\pmb{x_i}), \quad j=1,2,...,J \]
maka diperoleh \[ \boxed{ \pi_j(\pmb{x_i})= \frac{exp(\alpha_j+\pmb{x_i}^T\pmb{\beta})} {1+exp(\alpha_j+\pmb{x_i}^T\pmb{\beta})} - \frac{exp(\alpha_{j-1}+\pmb{x_i}^T\pmb{\beta})} {1+exp(\alpha_{j-1}+\pmb{x_i}^T\pmb{\beta})}, \quad j=1,2,...,J } \] dengan \(\frac{exp(\alpha_0+\pmb{x_i}^T\pmb{\beta})}{1+exp(\alpha_0+\pmb{x_i}^T\pmb{\beta})}=0\) dan \(\frac{exp(\alpha_{j-1}+\pmb{x_i}^T\pmb{\beta})}{1+exp(\alpha_{j-1}+\pmb{x_i}^T\pmb{\beta})}=0=1\).
Sebagai contoh dimisalkan variabel respon memiliki tiga buah kategori \((J=3)\), maka model regresi logistik ordinal yang terbentuk adalah \[ \begin{aligned} logit[P(Y_i\leq 1|\pmb{x_i})]&= ln \begin{bmatrix} \frac{P(Y_i\leq 1|\pmb{x_i})} {1-P(Y_i\leq 1|\pmb{x_i})} \end{bmatrix} = \alpha_1+\pmb{x_i}^T\pmb{\beta} \\ P(Y_i\leq 1|\pmb{x_i})&= \frac{exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})} {1+exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})} \\ logit[P(Y_i\leq 2|\pmb{x_i})]&= ln \begin{bmatrix} \frac{P(Y_i\leq 2|\pmb{x_i})} {2-P(Y_i\leq 2|\pmb{x_i})} \end{bmatrix} = \alpha_2+\pmb{x_i}^T\pmb{\beta} \\ P(Y_i\leq 2|\pmb{x_i})&= \frac{exp(\alpha_2+\pmb{x_i}^T\pmb{\beta})} {1+exp(\alpha_2+\pmb{x_i}^T\pmb{\beta})} \end{aligned} \]
Sehingga diperoleh peluang untuk masing-masing respon yaitu:
Peluang Kategori Pertama \[ \begin{aligned} \pi_1(\pmb{x_i})&=P(Y_i=1|\pmb{x_i})=P(Y_i\leq 1|\pmb{x_i}) \\ &= \frac{exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})} {1+exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})} \end{aligned} \]
Peluang Kategori Kedua \[ \begin{aligned} \pi_2(\pmb{x_i})&=P(Y_i=2|\pmb{x_i})=P(Y_i\leq 2|\pmb{x_i})-P(Y_i\leq 1|\pmb{x_i}) \\ &= \frac{exp(\alpha_2+\pmb{x_i}^T\pmb{\beta})} {1+exp(\alpha_2+\pmb{x_i}^T\pmb{\beta})}- \frac{exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})} {1+exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})} \end{aligned} \]
Peluang Kategori Ketiga \[ \begin{aligned} \pi_3(\pmb{x_i})&=P(Y_i=3|\pmb{x_i})=P(Y_i\leq 3|\pmb(x_i))-P(Y_i\leq 2|\pmb{x_i}) \\ &= 1-\frac{exp(\alpha_2+\pmb{x_i}^T\pmb{\beta})} {1+exp(\alpha_2+\pmb{x_i}^T\pmb{\beta})} \end{aligned} \]
Penaksiran parameter pada regresi logistik ordinal dilakukan dengan menggunakan metode Maximum Likelihood Estimation (MLE) (Agresti,2002). Jika diambil n sampel vektor variabel random \(\pmb{Y}_1, \pmb{Y}_2,..., \pmb{Y}_n\) dengan \(\pmb{Y}_i=\begin{bmatrix} y_{i1} & y_{i2} &...& y_{i,J-1}\end{bmatrix}^T\) yang berdistribusi multinomial dengan peluang hasil kategori ke-j adalah \(\pi_j(\pmb{x_i})\), maka terbentuk fungsi likelihood yaitu
\[ \begin{aligned} l(\pmb{\theta}) &= \prod_{i=1}^n \prod_{j=1}^J (\pi_j(\pmb{x_i}))^{y_{ij}} \\ &= \prod_{i=1}^n \prod_{j=1}^J \begin{bmatrix} \frac{exp(\alpha_j+\pmb{x_i}^T\pmb{\beta})} {1+exp(\alpha_j+\pmb{x_i}^T\pmb{\beta})} - \frac{exp(\alpha_{j-1}+\pmb{x_i}^T\pmb{\beta})} {1+exp(\alpha_{j-1}+\pmb{x_i}^T\pmb{\beta})} \end{bmatrix}^{y_{ij}} \end{aligned} \]
Menurut Hosmer dan Lemeshow (2000) prinsip dari metode MLE adalah mengestimasi vektor parameter \(\pmb{\theta}=\begin{bmatrix} \alpha_1 & \alpha_2 ... \alpha_{J-1} & \beta_1 & \beta_2 &... & \beta_{J-1} \end{bmatrix}^T\) dengan memaksimalkan fungsi likelihood lalu dilakukan transformasi ln pada fungsi likelihood sehingga terbentuk fungsi ln-likelihood sebagai berikut. \[ \begin{aligned} L(\pmb{\theta}) &= ln \begin{bmatrix} l(\pmb{\theta}) \end{bmatrix} \\ &= \prod_{i=1}^n \prod_{j=1}^J y_{ij} ln \begin{bmatrix} \frac{exp(\alpha_j+\pmb{x_i}^T\pmb{\beta})} {1+exp(\alpha_j+\pmb{x_i}^T\pmb{\beta})}- \frac{exp(\alpha_{j-1}+\pmb{x_i}^T\pmb{\beta})} {1+exp(\alpha_{j-1}+\pmb{x_i}^T\pmb{\beta})} \end{bmatrix} \end{aligned} \]
yang apabila dimisalkan terdapat tiga kategori pada variabel respon \((G=3)\), maka fungsi ln-likelihood menjadi \[ L(\pmb{\theta})=\sum_{i=1}^n \left\{y_{i1}ln \begin{bmatrix}\frac{exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})}{1+exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})}\end{bmatrix}+ y_{i2}ln\begin{bmatrix}\frac{exp(\alpha_2+\pmb{x_i}^T\pmb{\beta})}{1+exp(\alpha_2+\pmb{x_i}^T\pmb{\beta})}-\frac{exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})}{1+exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})}\end{bmatrix} + y_{i3}ln\begin{bmatrix}1-\frac{exp(\alpha_2+\pmb{x_i}^T\pmb{\beta})}{1+exp(\alpha_2+\pmb{x_i}^T\pmb{\beta})}\end{bmatrix} \right\} \]
yang dapat disederhanakan dalam bentuk \[ L(\pmb{\theta})=\sum_{i=1}^n \left\{ y_{i1}(\alpha_1+\pmb{x_i}^T\pmb{\beta})-(y_{i1}+y_{i2}ln\left[1+exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})\right]+y_{i2}ln\left[exp(\alpha_2+\pmb{x_i}^T\pmb{\beta})-exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})\right]+(y_{i1}-1)ln\left[1+exp(\alpha_2+\pmb{x_i}^T\pmb{\beta})\right] \right\} \]
Estimasi parameter dengan metode Maximum Likelihood Estimation (MLE) dilakukan dengan menurunkan fungsi log-likelihood terhadap parameter yang akan diestimasi. Turunan parsial pertama dari fungsi ln-likelihood tersebut kemudian disamakan dengan nol untuk memperoleh nilai parameter yang memaksimumkan fungsi likelihood. Oleh karena itu, langkah awal dalam proses estimasi adalah menentukan turunan parsial pertama dari fungsi ln-likelihood terhadap parameter yang diestimasi. $$ \[\begin{aligned} \frac{\partial L(\pmb{\theta})}{\partial \alpha_1}&= \sum_{i=1}^n \left\{ y_{i1}-(y_{i1}+y_{i2})\frac{exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})}{1+exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})}-y_{i2}\frac{exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})}{1+exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})-exp(\alpha_2+\pmb{x_i}^T\pmb{\beta})} \right\}=0 \\ \frac{\partial L(\pmb{\theta})}{\partial \alpha_2} &= \sum_{i=1}^n \left\{ y_{i2}\frac{exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})}{1+exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})-exp(\alpha_2+\pmb{x_i}^T\pmb{\beta})}+(y_{i1}-1)\frac{exp(\alpha_2+\pmb{x_i}^T\pmb{\beta})}{1+exp(\alpha_2+\pmb{x_i}^T\pmb{\beta})} \right\}=0 \\ \frac{\partial L(\pmb{\theta})}{\partial \beta}&= \sum_{i=1}^n \left\{ (y_{i1}+y_{i2})\pmb(x_i)^T\frac{1}{1+exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})}+(y_{i1}-1)\pmb{x_i}^T\frac{exp(\alpha_2+\pmb{x_i}^T\pmb{\beta})}{1+exp(\alpha_2+\pmb{x_i}^T\pmb{\beta})} \right\}=0 \end{aligned}\]$$
Turunan parsial pertama dari fungsi ln-likelihood yang akan diestimasi umumnya menghasilkan persamaan yang bersifat nonlinear terhadap parameter. Oleh karena itu, penentuan estimasi parameter dengan metode kemungkinan maksimum tidak dapat diselesaikan secara langsung dengan cara analitik sederhana. Proses estimasi biasanya memerlukan prosedur numerik yang bersifat iteratif untuk memperoleh nilai parameter yang konvergen. Salah satu metode iteratif yang sering digunakan adalah metode Newton–Raphson, yaitu algoritma numerik yang mencari solusi persamaan nonlinear dengan mendekati titik dimana fungsi mencapai nilai maksimum atau minimum. Metode ini bekerja dengan memperbarui nilai parameter secara berulang berdasarkan turunan pertama dan kedua dari fungsi yang dioptimalkan sampai diperoleh nilai yang stabil.
Oleh karena itu diperlukan turunan parsial kedua fungsi ln-likelihood terhadap parameter yang akan diestimasi yaitu sebagai berikut. \[ \begin{aligned} \frac{\partial^2 L(\pmb{\theta})}{\partial \alpha_1^2}&= \sum_{i=1}^n \left\{ -(y_{i1}+y_{i2})\frac{exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})}{\left[1+exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})\right]^2}-y_{i2}\frac{exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})exp(\alpha_2+\pmb{x_i}^T\pmb{\beta})}{\left[exp(\alpha_2+\pmb{x_i}^T\pmb{\beta})-exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})\right]^2} \right\} \\ \frac{\partial^2 L(\pmb{\theta})}{\partial \alpha_1\partial\alpha_2} &= \sum_{i=1}^n y_{i2}\frac{exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})exp(\alpha_2+\pmb{x_i}^T\pmb{\beta})}{\left[exp(\alpha_2+\pmb{x_i}^T\pmb{\beta})-exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})\right]^2} \\ \frac{\partial^2 L(\pmb{\theta})}{\partial \alpha_1\partial\beta} &= \sum_{i=1}^n \left\{ -(y_{i1}+y_{i2})\pmb{x_i}^T\frac{exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})}{\left[1+exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})\right]^2} \right\} \\ \frac{\partial^2 L(\pmb{\theta})}{\partial \alpha_2^2}&= \sum_{i=1}^n \left\{ -y_{i2}\frac{exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})exp(\alpha_2+\pmb{x_i}^T\pmb{\beta})}{\left[exp(\alpha_2+\pmb{x_i}^T\pmb{\beta}-exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})\right]^2}+(y_{i1}-1)\frac{exp(\alpha_2+\pmb{x_i}^T\pmb{\beta})}{\left[1+exp(\alpha_2+\pmb{x_i}^T\pmb{\beta})\right]^2} \right\} \\ \frac{\partial^2 L(\pmb{\theta})}{\partial \alpha_2\partial\beta} &= \sum_{i=1}^n \left\{ (y_{i1}-1)\pmb{x_i}^T\frac{exp(\alpha_2+\pmb{x_i}^T\pmb{\beta})}{\left[1+exp(\alpha_2+\pmb{x_i}^T\pmb{\beta})\right]^2} \right\} \\ \frac{\partial^2 L(\pmb{\theta})}{\partial \beta \partial{\beta}^T}&= \sum_{i=1}^n \left\{ -(y_{i1}+y_{i2})\pmb{x_i}^T\pmb{x_i}\frac{exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})}{\left[1+exp(\alpha_1+\pmb{x_i}^T\pmb{\beta})\right]^2}+(y_{i1}-1)\pmb{x_i}^T\pmb{x_i}\frac{exp(\alpha_2+\pmb{x_i}^T\pmb{\beta})}{\left[1+exp(\alpha_2+\pmb{x_i}^T\pmb{\beta})\right]^2} \right\} \end{aligned} \]
Menurut Alan Agresti (2002), persamaan yang digunakan dalam proses iterasi Newton–Raphson untuk memperoleh nilai estimasi parameter \(\pmb{\theta}\) dinyatakan sebagai \[ \pmb{\theta}^{(t+1)})=\pmb{\theta}^{(t)}-\left[\pmb{H}(\pmb{\theta}^{(t)})\right]^{-1}\pmb{q}(\pmb{\theta}^{(t)}) \]
dengan \(\pmb{H}(\pmb{\theta})\) merupakan matriks nonsingular yang elemen-elemennya berupa turunan parsial kedua dari fungsi ln-likelihood terhadap parameter yang diestimasi. Sementara itu, \(\pmb{q}(\pmb{\theta})\) adalah vektor yang berisi turunan parsial pertama dari fungsi log-likelihood terhadap parameter yang diestimasi. Notasi t menyatakan banyaknya iterasi (\(t=0,1,2,...\)). Dengan demikian, bentuk elemen dari vektor \(\pmb{q}(\pmb{\theta})\) dan matriks \(\pmb{H}(\pmb{\theta})\) dapat dituliskan sebagai berikut. \[ \begin{aligned} \pmb{q}(\pmb{\theta}) &= \begin{bmatrix} \frac{\partial L(\pmb{\theta})}{\partial \alpha_1} & \frac{\partial L(\pmb{\theta})}{\partial \alpha_2} & \frac{\partial L(\pmb{\theta})}{\partial \pmb{\beta}} \end{bmatrix} \\ \pmb{H}(\pmb{\theta}) &= \begin{bmatrix} \frac{\partial^2 L(\pmb{\theta})}{\partial \alpha_1^2} & \frac{\partial^2 L(\pmb{\theta})}{\partial \alpha_1 \partial \alpha_2} & \frac{\partial^2 L(\pmb{\theta})}{\partial \alpha_1 \partial \pmb{\beta}} \\ \frac{\partial^2 L(\pmb{\theta})}{\partial \alpha_1 \partial \alpha_2} & \frac{\partial^2 L(\pmb{\theta})}{\partial \alpha_2^2} & \frac{\partial^2 L(\pmb{\theta})}{\partial \alpha_2 \partial \beta} \\ \frac{\partial^2 L(\pmb{\theta})}{\partial \alpha_1 \partial \pmb{\beta}} & \frac{\partial^2 L(\pmb{\theta})}{\partial \alpha_2 \partial \pmb{\beta}} & \frac{\partial^2 L(\pmb{\theta})}{\partial \pmb{\beta} \partial \pmb{\beta}}^T\\ \end{bmatrix} \end{aligned} \]
Pada regresi logistik ordinal terdapat asumsi penting yaitu asumsi proportional odds atau parallel lines. Asumsi ini menyatakan bahwa pengaruh variabel independen terhadap logit kumulatif bersifat konstan pada setiap batas kategori variabel respon. Dengan kata lain, koefisien regresi yang dihasilkan untuk setiap kategori respon adalah sama.
Pengujian asumsi ini dapat dilakukan menggunakan Brant Test. Hipotesis yang digunakan dalam pengujian ini adalah sebagai berikut:
\(H_0\) : Asumsi proportional odds terpenuhi.
\(H_1\) : Asumsi proportional odds tidak terpenuhi.
Kriteria pengambilan keputusan dilakukan dengan membandingkan nilai p-value dengan taraf signifikansi 0,05. Jika p-value lebih besar dari 0,05 maka \(H_0\) tidak ditolak sehingga dapat disimpulkan bahwa asumsi proportional odds terpenuhi.
Dalam perangkat lunak R, pengujian ini dapat dilakukan menggunakan
Brant Test yang tersedia pada paket brant Setelah
model regresi logistik ordinal diestimasi menggunakan fungsi
polr() dari paket MASS.
Uji Goodness of Fit dilakukan untuk mengetahui apakah model regresi logistik ordinal yang dibangun telah sesuai dengan data yang digunakan. Salah satu metode yang digunakan adalah uji Deviance yang membandingkan model yang diestimasi dengan model jenuh.
Hipotesis:
\(H_0\) : Model sesuai dengan data (good fit).
\(H_1\) : Model tidak sesuai dengan data.
Statistik Uji Deviance \[ D=-2\sum_{i=1}^n \sum_{j=1}^J \left[ y_{ij}ln \left( \frac{\hat{\pi_{ij}}}{y_{ij}} \right) \right] \]
dengan \(\hat{\pi_{ij}}=\hat{\pi_g}(x_i)\) merupakan peluang observasi ke-i pada kategori ke-j
Apabila nilai \(Deviance >
\chi^2_{(J-(p+1);\alpha)}\) atau nilai \(p-value < 0,05\) maka tidak ada cukup
bukti untuk menerima \(H_0\) sehingga
model dianggap tidak memiliki kecocokan yang baik terhadap data. Nilai
Deviance juga dapat diperoleh dari output pada fungsi
polr() pada perangkat lunak R.
Uji Rasio Likelihood diperoleh dengan cara membandingkan fungsi ln-likelihood dari seluruh variabel bebas dengan fungsi ln-likelihood tanpa variabel bebas. Uji Rasio Likelihood digunakan untuk mengetahui apakah variabel bebas yang terdapat dalam model berpengaruh nyata atau tidak secara keseluruhan.
Hipotesis
\(H_0:\beta_1=\beta_2=...=\beta_p\)
\(H_1\): Paling sedikit salah satu dari \(\beta_k\neq 0\) dengan \(k=1,2,...,p\)
Statistik Uji Rasio Likelihood \[ G^2 =-2ln\left(\frac{\text{likelihood tanpa variabel bebas}}{\text{likelihood dengan variabel bebas}}\right) \]
Kriteria uji pada uji adalah apabila \(H_0\) ditolak jika nilai \(G^2 > \chi^2_{(\alpha,p)}\) atau nilai
\(p-value < \alpha\). Nilai
statistik hitung \(G^2\) juga dapat
diperoleh menggunakan fungsi pR2 pada package
pscl.
Uji Wald digunakan untuk mengetahui apakah masing-masing variabel independen berpengaruh secara signifikan terhadap variabel dependen dalam model regresi logistik ordinal. Uji ini dilakukan secara parsial untuk setiap koefisien parameter yang terdapat dalam model.
Hipotesis yang digunakan dalam uji Wald adalah: \(H_0:\beta_k=0\)
\(H_1:\beta_k\neq 0\) dengan \(k=1,2,...,p\)
Statistik Uji Rasio Likelihood \[ W_k= \left[ \frac{\hat{\beta_k}}{SE(\hat{\beta_k})} \right]^2 \]
Kita dapat menolak \(H_0\) apabila nilai \(W_k>\chi^2_{(\alpha;1)}\) atau nilai \(p-value<\alpha(0,05)\) sehingga kesimpulan yang diperoleh adalah variabel independen tersebut berpengaruh signifikan terhadap variabel dependen.
Dalam perangkat lunak R, uji Wald dapat diperoleh langsung dari hasil
estimasi model regresi logistik ordinal menggunakan fungsi
polr() pada paket MASS. Nilai statistik uji dapat dilihat
pada output fungsi summary() yang menampilkan nilai
t-value untuk setiap koefisien parameter.
Data yang digunakan dalam studi kasus ini adalah dataset
wine yang tersedia pada paket ordinal di perangkat
lunak R. Dataset ini sering digunakan sebagai contoh dalam analisis
regresi logistik ordinal karena variabel responnya berbentuk kategori
berurutan (ordinal).
Dataset wine berisi hasil penilaian kualitas anggur oleh sejumlah panelis terhadap anggur yang diproses dengan perlakuan tertentu. Variabel respon dalam dataset ini adalah rating, yang menunjukkan tingkat penilaian kualitas anggur oleh panelis. Variabel ini bersifat ordinal karena terdiri dari beberapa kategori penilaian yang memiliki urutan tingkat kualitas, mulai dari kualitas rendah hingga kualitas yang lebih tinggi.
Dalam penelitian ini digunakan dua variabel independen, yaitu temp dan contact. Variabel temp menunjukkan suhu proses fermentasi anggur, sedangkan variabel contact menunjukkan lama kontak antara kulit anggur dengan sari buah selama proses pembuatan anggur. Kedua variabel tersebut merupakan faktor proses produksi yang secara teoritis dapat memengaruhi karakteristik dan kualitas akhir anggur.
Pemilihan variabel temp dan contact sebagai variabel independen dalam model dilakukan karena keduanya merupakan faktor yang secara langsung berkaitan dengan proses pembuatan anggur. Perbedaan suhu fermentasi dapat memengaruhi reaksi kimia selama proses produksi, sementara lama kontak dengan kulit anggur dapat memengaruhi warna, aroma, dan cita rasa anggur. Oleh karena itu, variasi pada kedua faktor tersebut diduga dapat memengaruhi penilaian kualitas anggur yang diberikan oleh panelis.
library(MASS)
library(dplyr)
library(ordinal)
library(kableExtra)
# Input Data
ordinaldf<- ordinal::wine[,2:4]Berikut adalah preview dataset yang digunakan pada studi kasus ini.
| rating | temp | contact |
|---|---|---|
| 2 | cold | no |
| 3 | cold | no |
| 3 | cold | yes |
| 4 | cold | yes |
| 4 | warm | no |
| 4 | warm | no |
| 5 | warm | yes |
| 5 | warm | yes |
| 1 | cold | no |
| 2 | cold | no |
Sebelum melakukan analisis regresi logistik ordinal, beberapa
variabel pada dataset perlu disesuaikan tipe datanya agar sesuai dengan
kebutuhan model. Salah satu cara yang dilakukan adalah dengan
menggunakan fungsi as.factor() di R.
Fungsi as.factor() digunakan untuk mengubah tipe data
suatu variabel menjadi faktor, yaitu tipe data kategorik yang terdiri
dari beberapa kategori atau level. Transformasi ini penting dilakukan
apabila variabel yang digunakan dalam analisis merupakan variabel
kategorik tetapi masih tersimpan dalam bentuk numerik atau karakter.
Model regresi logistik ordinal pada penelitian ini dibentuk
menggunakan fungsi polr() dari paket MASS
dalam perangkat lunak R. Model tersebut dituliskan sebagai berikut:
Perintah tersebut digunakan untuk mengestimasi model regresi logistik
ordinal dengan variabel rating sebagai variabel dependen dan
variabel temp serta contact sebagai variabel
independen. Argumen method = "logistic" menunjukkan bahwa
model yang digunakan adalah model regresi logistik ordinal dengan fungsi
tautan logit.
Berdasarkan hasil estimasi model yang diperoleh, maka model peluang kumulatif untuk setiap kategori rating dapat dituliskan sebagai berikut.
## Call:
## polr(formula = rating ~ temp + contact, data = ordinaldf, method = "logistic")
##
## Coefficients:
## Value Std. Error t value
## tempwarm 2.503 0.5287 4.735
## contactyes 1.528 0.4766 3.205
##
## Intercepts:
## Value Std. Error t value
## 1|2 -1.3444 0.5171 -2.5998
## 2|3 1.2508 0.4379 2.8565
## 3|4 3.4669 0.5978 5.7998
## 4|5 5.0064 0.7309 6.8496
##
## Residual Deviance: 172.9838
## AIC: 184.9838Model untuk \(P(Y_i\leq 1|x_i)\) \[ \begin{aligned} P(Y_i\leq 1|x_i) &= \frac{e^{\alpha_1+\beta_1(temp_warm)+\beta_2(contact_yes)}}{1+e^{\alpha_1+\beta_1(temp_warm)+\beta_2(contact_yes)}} \\ P(Y_i\leq 1|x_i) &= \frac{e^{-1.344+2.503(temp_warm)+1.528(contact_yes)}}{1+e^{-1.344+2.503(temp_warm)+1.528(contact_yes)}} \end{aligned} \]
Model untuk \(P(Y_i\leq 2|x_i)\) \[ \begin{aligned} P(Y_i\leq 2|x_i) &= \frac{e^{\alpha_2+\beta_1(temp_warm)+\beta_2(contact_yes)}}{1+e^{\alpha_2+\beta_1(temp_warm)+\beta_2(contact_yes)}} \\ P(Y_i\leq 2|x_i) &= \frac{e^{1.2508+2.503(temp_warm)+1.528(contact_yes)}}{1+e^{1.2508+2.503(temp_warm)+1.528(contact_yes)}} \end{aligned} \]
Model untuk \(P(Y_i\leq 3|x_i)\) \[ \begin{aligned} P(Y_i\leq 3|x_i) &= \frac{e^{\alpha_3+\beta_1(temp_warm)+\beta_2(contact_yes)}}{1+e^{\alpha_3+\beta_1(temp_warm)+\beta_2(contact_yes)}} \\ P(Y_i\leq 3|x_i) &= \frac{e^{3.467+2.503(temp_warm)+1.528(contact_yes)}}{1+e^{3.467+2.503(temp_warm)+1.528(contact_yes)}} \end{aligned} \]
Model untuk \(P(Y_i\leq 4|x_i)\) \[ \begin{aligned} P(Y_i\leq 4|x_i) &= \frac{e^{\alpha_4+\beta_1(temp_warm)+\beta_2(contact_yes)}}{1+e^{\alpha_4+\beta_1(temp_warm)+\beta_2(contact_yes)}} \\ P(Y_i\leq 4|x_i) &= \frac{e^{5.0064+2.503(temp_warm)+1.528(contact_yes)}}{1+e^{5.0064+2.503(temp_warm)+1.528(contact_yes)}} \end{aligned} \]
Untuk mengetahui apakah pengaruh variabel independen terhadap logit kumulatif bersifat konstan pada setiap batas kategori variabel respon perlu dilakukan pengujian asumsi Proportional Odds yang dilakukan menggunakan Brant Test seperti di bawah ini.
## --------------------------------------------
## Test for X2 df probability
## --------------------------------------------
## Omnibus 1.52 6 0.96
## tempwarm 0.93 3 0.82
## contactyes 0.86 3 0.84
## --------------------------------------------
##
## H0: Parallel Regression Assumption holdsBerdasarkan hasil uji Brant yang diperoleh, nilai p-value pada uji omnibus sebesar 0,96. Karena nilai tersebut lebih besar dari taraf signifikansi yang digunakan (\(\alpha = 0,05\)), maka \(H_0\) gagal ditolak. Hal ini menunjukkan bahwa secara keseluruhan model regresi logistik ordinal yang digunakan telah memenuhi asumsi proportional odds.
Selain itu, pengujian juga dilakukan pada masing-masing variabel
independen. Variabel temp_warm memiliki p-value
sebesar 0,82, sedangkan variabel contac_tyes memiliki
p-value sebesar 0,84. Kedua nilai tersebut juga lebih besar
dari 0,05, sehingga dapat disimpulkan bahwa tidak terdapat pelanggaran
asumsi proportional odds pada masing-masing variabel
independen.
Dengan demikian, dapat disimpulkan bahwa model regresi logistik ordinal yang dibentuk telah memenuhi asumsi parallel regression, sehingga model tersebut layak digunakan untuk analisis lebih lanjut.
Uji kecocokan model (goodness of fit) dilakukan untuk mengetahui apakah model regresi logistik ordinal yang dibentuk telah sesuai dengan data yang digunakan.
## [1] 172.9838Diperoleh Deviance adalah sebesar 172.9839 dengan nilai \(\chi^2_{(\alpha,J-(p+1))}\) sebesar 5.9914645. Sehingga diperoleh kesimpulan jika \(H_0\) ditolak karena nilai \(Deviance (172.9838467) > \chi^2_{(\alpha,J-(p+1))}(5.9914645)\) yang berarti model sesuai atau tidak ada perbedaan antara observasi dan prediksi.
Berdasarkan hasil perhitungan menggunakan fungsi pR2(),
diperoleh beberapa ukuran Pseudo R² yang digunakan untuk
menilai kemampuan model regresi logistik ordinal dalam menjelaskan
variabel respon.
## fitting null model for pseudo-r2## llh llhNull G2 McFadden r2ML r2CU
## -86.4919234 -103.7190774 34.4543080 0.1660944 0.3803076 0.4028997Statistik Uji: \[ \begin{aligned} G^2 &=-2ln\left(\frac{\text{likelihood tanpa variabel bebas}}{\text{likelihood dengan variabel bebas}}\right) \\ G^2 &= -2ln \left(\frac{llhNull}{llh}\right) \\ G^2 &= -2ln(llhNull)-(-2ln(llh)) \\ G^2 &= -2(-103.7190774)+2(-86.4919234) \\ G^2 &= 34.454308 \end{aligned} \]
dengan nilai \(\chi^2_{(\alpha,p)}\) dengan \(\alpha=0,05\) dan p= banyaknya variabel bebas = 2 (temp & contact) diperoleh nilai sebesar 5.9914645. Didapatkan jika nilai statistik hitung \(G^2 (34.454308)>\chi^2_{(\alpha,p)}(5.9914645)\) sehingga \(H_0\) ditolak sehingga dapat disimpulkan jika secara bersama-sama variabel bebas temp & contact memengaruhi model.
Uji Wald digunakan untuk mengetahui apakah masing-masing variabel independen dalam model regresi logistik ordinal berpengaruh secara signifikan terhadap variabel dependen. Pengujian ini dilakukan secara parsial terhadap setiap koefisien regresi yang diestimasi dalam model.
coef_estimates<- coef(model_polr)
std_errors<- sqrt(diag(vcov(model_polr)))[names(coef(model_polr))]
z_scores<- coef_estimates/std_errors
p_value<- 2* (1-pnorm(abs(z_scores)))
wald_test_results<- data.frame(Coefficients=coef_estimates,
Std_Errors=std_errors,
Z_Scores=z_scores,
P_values=p_value)
print(wald_test_results)## Coefficients Std_Errors Z_Scores P_values
## tempwarm 2.503073 0.5286779 4.734589 2.194990e-06
## contactyes 1.527786 0.4766213 3.205450 1.348517e-03coef_estimates<- coef(model_polr)
std_errors<- sqrt(diag(vcov(model_polr)))[names(coef(model_polr))]
z_scores<- coef_estimates/std_errors
p_value<- 2* (1-pnorm(abs(z_scores)))
wald_test_results<- data.frame(Coefficients=coef_estimates,
Std_Errors=std_errors,
Z_Scores=z_scores,
P_values=p_value)
print(wald_test_results)## Coefficients Std_Errors Z_Scores P_values
## tempwarm 2.503073 0.5286779 4.734589 2.194990e-06
## contactyes 1.527786 0.4766213 3.205450 1.348517e-03Variabel temp_warm memiliki nilai koefisien sebesar
2.503073, nilai Z-score sebesar 4.734589, dan p-value
sebesar 2.194990 × 10⁻⁶. Karena nilai p-value lebih kecil dari
taraf signifikansi 0,05, maka H₀ ditolak. Hal ini menunjukkan bahwa
variabel temp berpengaruh signifikan terhadap rating. Koefisien
yang bernilai positif menunjukkan bahwa penggunaan suhu fermentasi warm
meningkatkan peluang anggur memperoleh kategori rating yang lebih
tinggi.
Selanjutnya, variabel contact_yes memiliki nilai
koefisien sebesar 1.527786, nilai Z-score sebesar 3.205450, dan
p-value sebesar 0.001348517. Karena nilai p-value juga
lebih kecil dari 0,05, maka H₀ ditolak. Hal ini menunjukkan bahwa
variabel contact berpengaruh signifikan terhadap rating. Nilai
koefisien yang positif menunjukkan bahwa adanya kontak antara kulit
anggur dan sari buah selama proses produksi meningkatkan peluang anggur
memperoleh rating yang lebih tinggi.
Dengan demikian, berdasarkan hasil uji Wald dapat disimpulkan bahwa variabel temp dan contact secara parsial berpengaruh signifikan terhadap rating dalam model regresi logistik ordinal yang dibentuk.
Odds ratio digunakan untuk mengetahui besarnya perubahan peluang suatu observasi berada pada kategori respon yang lebih tinggi akibat perubahan pada variabel independen. Nilai odds ratio diperoleh dengan melakukan eksponensial terhadap koefisien regresi.
## OR 2.5 % 97.5 %
## tempwarm 12.219985 4.525654 36.42340
## contactyes 4.607962 1.851150 12.09059Berdasarkan hasil perhitungan, variabel temp_warm
memiliki nilai odds ratio sebesar 12.219985 dengan interval kepercayaan
95% sebesar 4.525654 hingga 36.42340. Hal ini menunjukkan bahwa anggur
yang diproses dengan suhu fermentasi warm memiliki peluang
sekitar 12,22 kali lebih besar untuk memperoleh rating yang
lebih tinggi dibandingkan dengan suhu referensi (cold). Selain
itu, karena interval kepercayaan tidak mencakup nilai 1, maka dapat
disimpulkan bahwa pengaruh variabel temp terhadap rating bersifat
signifikan.
Selanjutnya, variabel contact_yes memiliki nilai odds
ratio sebesar 4.607962 dengan interval kepercayaan 95% sebesar 1.851150
hingga 12.09059. Hal ini menunjukkan bahwa anggur yang mengalami kontak
dengan kulit anggur memiliki peluang sekitar 4,61 kali lebih besar untuk
memperoleh rating yang lebih tinggi dibandingkan dengan anggur
yang tidak mengalami kontak. Interval kepercayaan yang juga tidak
mencakup nilai 1 menunjukkan bahwa variabel contact memiliki pengaruh
yang signifikan terhadap rating.
Berdasarkan hasil analisis regresi logistik ordinal yang telah dilakukan terhadap dataset wine, dapat disimpulkan bahwa variabel temp dan contact berpengaruh terhadap rating anggur yang diberikan oleh panelis.
Hasil pengujian asumsi menggunakan Brant Test menunjukkan bahwa model memenuhi asumsi proportional odds, sehingga model regresi logistik ordinal layak digunakan dalam analisis. Selanjutnya, hasil Likelihood Ratio Test menunjukkan bahwa model yang memasukkan variabel temp dan contact secara signifikan lebih baik dibandingkan dengan model tanpa variabel prediktor.
Berdasarkan uji Wald, kedua variabel independen yaitu temp dan contact terbukti berpengaruh signifikan terhadap rating. Selain itu, hasil perhitungan odds ratio menunjukkan bahwa penggunaan suhu fermentasi warm meningkatkan peluang anggur memperoleh kategori rating yang lebih tinggi sekitar 12,22 kali dibandingkan suhu referensi (cold). Sementara itu, adanya kontak antara kulit anggur dengan sari buah meningkatkan peluang memperoleh rating yang lebih tinggi sekitar 4,61 kali dibandingkan tanpa kontak.
Nilai Pseudo R² yang diperoleh menunjukkan bahwa model memiliki kemampuan yang cukup baik dalam menjelaskan variasi pada variabel rating. Dengan demikian, dapat disimpulkan bahwa faktor suhu fermentasi dan kontak kulit anggur dalam proses produksi berperan dalam memengaruhi kualitas anggur yang dinilai oleh panelis.
How to Fit an Ordinal Regression Model Using R
Ordinal Logistic Regression|R Data Analysis Examples
[Wulandari, Yuciana. Regresi Logistik Ordinal]