Bagaimana agar sesuai dengan regresi seperti di R?

9

Saya punya beberapa data deret waktu di mana variabel yang diukur adalah bilangan bulat positif (jumlah). Saya ingin menguji apakah ada tren naik dari waktu ke waktu (atau tidak). Variabel independen (x) berada dalam kisaran 0-500 dan variabel dependen (y) berada dalam kisaran 0-8.

Saya pikir saya menjawab ini dengan menyesuaikan regresi formulir y = floor(a*x + b)menggunakan ordinary least square (OLS).

Bagaimana cara saya melakukan ini menggunakan R (atau Python)? Apakah ada paket yang ada untuknya, atau lebih baik saya menulis algoritma saya sendiri?

PS: Saya tahu ini bukan teknik yang ideal, tetapi saya perlu melakukan analisis yang relatif sederhana yang sebenarnya bisa saya pahami - latar belakang saya adalah biologi, bukan matematika. Saya tahu saya melanggar asumsi tentang kesalahan dalam variabel yang diukur, dan independensi pengukuran dari waktu ke waktu.

r regression python afaulconbridge
sumber

5

Meskipun secara matematis wajar untuk mencoba regresi dari bentuk ini, di belakangnya bersembunyi kesalahan statistik: istilah kesalahan sekarang akan sangat berkorelasi dengan nilai yang diprediksi. Itu pelanggaran asumsi OLS yang cukup kuat. Sebagai gantinya, gunakan teknik berbasis hitungan seperti yang disarankan oleh balasan Greg Snow. (Namun, saya dengan senang hati mengangkat pertanyaan ini, karena itu mencerminkan beberapa pemikiran dan kepintaran nyata. Terima kasih telah mengajukannya di sini!)

whuber

11

Anda dapat memasukkan model yang Anda nyatakan menggunakan fungsi nls(non-linear least square) R, tetapi seperti yang Anda katakan akan melanggar banyak asumsi dan masih mungkin tidak masuk akal (Anda mengatakan hasil yang diprediksi adalah acak di sekitar langkah fungsi, bukan nilai integer di sekitar hubungan yang meningkat dengan lancar).

Cara yang lebih umum untuk mencocokkan data hitung menggunakan regresi Poisson menggunakan glmfungsi dalam R, contoh pertama pada halaman bantuan adalah regresi Poisson, meskipun jika Anda tidak terbiasa dengan statistik, sebaiknya berkonsultasi dengan ahli statistik untuk memastikan Anda melakukan sesuatu dengan benar.

Jika nilai 8 adalah maksimum absolut (tidak mungkin untuk melihat jumlah yang lebih tinggi, bukan hanya itu yang Anda lihat) maka Anda dapat mempertimbangkan regresi logistik odds proporsional, ada beberapa alat untuk melakukan ini dalam paket R, tetapi Anda benar-benar harus melibatkan ahli statistik jika Anda ingin melakukan ini.

Greg Snow
sumber

"Anda mengatakan hasil yang diprediksi adalah acak di sekitar fungsi langkah, bukan nilai integer di sekitar hubungan yang meningkat dengan lancar" --- Itu adalah sesuatu yang belum saya pertimbangkan. Pada akhirnya, saya menggunakan regresi Poisson oleh glm. Ini bukan pilihan yang sempurna, tetapi "cukup baik" untuk apa yang saya butuhkan.

afaulconbridge

10

$\def\lf{\lfloor}\def\rf{\rfloor}\def\pnorm{\mathrm{pnorm}}$ Jelas bahwa saran Greg adalah hal pertama yang harus dicoba: Regresi Poisson adalah model alami dalam banyak beton. situasi.

Namun model yang Anda sarankan dapat terjadi misalnya ketika Anda mengamati data bulat: dengan kesalahan normal iid .

Y_{i} = ⌊ a x_{i} + b + ϵ_{i} ⌋,

$Y_i = \lf ax_i + b + \epsilon_i \rf,$

ϵ_{i}

$\epsilon_i$

Saya pikir ini menarik untuk melihat apa yang bisa dilakukan dengannya. Saya tunjukkan dengan cdf dari variabel normal standar. Jika , maka menggunakan notasi komputer yang dikenal. $F$ $\epsilon \sim \mathcal N(0,\sigma^2)$

\begin{aligned} P (⌊ a x + b + ϵ ⌋ = k) & = F (\frac{k - b + 1 - a x}{σ}) - F (\frac{k - b - a x}{σ}) \\ = p n o r m (k + 1 - a x - b, s d = σ) - p n o r m (k - a x - b, s d = σ), \end{aligned}

$\begin{align*} \mathbb P\left(\lf ax + b + \epsilon \rf = k\right) &= F\left({k-b+1-ax\over \sigma}\right) - F\left({k-b-ax\over \sigma}\right)\\ &= \pnorm(k+1-ax-b,sd=\sigma) - \pnorm(k-ax-b,sd=\sigma),\end{align*}$

Anda mengamati titik data . Kemungkinan log diberikan oleh Ini tidak identik dengan kuadrat terkecil. Anda dapat mencoba memaksimalkan ini dengan metode numerik. Berikut ini adalah ilustrasi dalam R: $(x_i,y_i)$

ℓ (a, b, σ) = \sum_{i} \log (F (\frac{y_{i} - b + 1 - a x_{i}}{σ}) - F (\frac{y_{i} - b - a x_{i}}{σ})) .

$\ell(a,b,\sigma) = \sum_i \log\left( F\left({y_i-b+1-ax_i\over \sigma}\right) - F\left({y_i-b-ax_i\over \sigma}\right) \right).$

log_lik <- function(a,b,s,x,y)
  sum(log(pnorm(y+1-a*x-b, sd=s) - pnorm(y-a*x-b, sd=s)));

x <- 0:20
y <- floor(x+3+rnorm(length(x), sd=3))
plot(x,y, pch=19)
optim(c(1,1,1), function(p) -log_lik(p[1], p[2], p[3], x, y)) -> r
abline(r$par[2], r$par[1], lty=2, col="red")
t <- seq(0,20,by=0.01)
lines(t, floor( r$par[1]*t+r$par[2]), col="green")

lm(y~x) -> r1
abline(r1, lty=2, col="blue");

model linear bulat

Dalam merah dan biru, garis ditemukan oleh maksimalisasi numerik dari kemungkinan ini, dan kuadrat terkecil, masing-masing. Tangga hijau adalah untuk ditemukan dari kemungkinan maksimum ... ini menunjukkan bahwa Anda dapat menggunakan kuadrat terkecil, hingga terjemahan oleh 0,5, dan mendapatkan hasil yang kira-kira sama; atau, bahwa kuadrat terkecil cocok dengan model mana adalah bilangan bulat terdekat. Data bulat begitu sering bertemu sehingga saya yakin ini diketahui dan telah dipelajari secara luas ... $ax+b$ $\lf ax +b\rf$ $a,b$ $b$

Y_{i} = [a x_{i} + b + ϵ_{i}],

$Y_i = [ a x_i + b +\epsilon_i],$

[x] = ⌊ x + 0.5 ⌋

$[x] = \lf x + 0.5 \rf$

Elvis
sumber

4

+1 Saya suka teknik ini dan benar-benar mengirimkan makalah ke jurnal analisis risiko beberapa tahun yang lalu. (Beberapa analis risiko cukup tertarik pada data bernilai interval.) Itu ditolak karena "terlalu matematis" untuk audiens mereka. :-(. Satu kiat: saat menggunakan metode numerik, selalu merupakan ide bagus untuk memberikan nilai awal yang baik untuk solusi. Pertimbangkan menerapkan OLS ke data mentah untuk mendapatkan nilai-nilai itu, lalu "poles" dengan pengoptimal angka.

whuber

Ya, ini saran yang bagus. Sebenarnya, dalam hal ini saya memilih nilai jarak jauh untuk menekankan bahwa "itu bekerja", tetapi dalam praktiknya saran Anda akan menjadi satu-satunya solusi untuk menghindari mulai dari daerah yang sangat datar, tergantung pada data ...

Elvis

Bagaimana agar sesuai dengan regresi seperti di R?

Jawaban: