Zero-truncated Poisson regression Опции

[ИНО]

Печальный у Вас жизненный опыт, но не стоит обобщать. По ссылкам ходят и читают, если источник доступен. Если недоступен, ходят, чтобы убедиться в этом. Но первое необходимое условие для этого - наличие этих-самых ссылок в удобоваримом виде. Интересно, сколько есть в природе книг, написанных некоим Ивановым Джонсоном с соавторами с количеством страниц более 188? Прикажете все проверить? И такие с позволения сказать ссылки я встречаю в Ваших постах уже далеко не первый раз.

Сообщение отредактировал ИНО - 10.11.2022 - 07:32

[100$]

Вопрос: а почему с аналогичной целью нельзя применить обычную регрессию Пуассона, просто предварительно вычтя из всех значений зависимой переменной единицу? Пусть такая модель предсказывает не общее количество счетных единиц, а количество единиц, добавленных к одной обязательно имеющейся. Потом, если надо, добавляем к предсказанной величине единичку, и дело в шляпе! Но, судя по тому, как извращается народ именно с Zero-truncated Poisson regression, и ищет ее программные реализации, очевидно, предложенный мною альтернативный подход неправомерен. Но сколько я не ломал голову, так и не понял, почему. Прошу более головастых подсказать.

Надо думать патамушта моделирование суммы X1+X2+,...,+Xn - An независимых одинаково распределенных величин требует перехода к т.н. обобщенному распределению Пуассона (ОПР).
Число n слагаемых в этой сумме является Пуассоновой величиной с параметром lambda.
Как любое безгранично делимое распределение (а ОПР является таковым) каждое обобщенное распределение Пуассона является пределом других ОПР, в т.ч. и "сдвинутых" на некоторую неслучайную алгебраическую величину An. (т.е. характеристическая ф-ция ОПР допускает такой "сдвиг").

Разумеется, проверка сказанного с помощью "обычной регрессии Пуассона" приведет к появлению жалобы типа

Я пробовал - для того же набора данных данных получаются существенно другие p, значит, очевидно, есть какой-то теоретический запрет для этого моего ноу-хау.
Понять бы, какой.

с единичной вероятностью.

Сообщение отредактировал 100$ - 12.11.2022 - 18:14

[ИНО]

Благодарю за первый ответ качественный по существу проблемы! Но непонятно: почему мы должны подходить к величине X1 (т. е. первому прыщу из описанного примера) как к случайной, а не как к постоянной? Ведь этот первый прыщ есть примененное условие включения пациента в выборку. Не вижу причин не принимать X1 константой, равной 1, и моделировать только X2.

[100$]

...Но непонятно: почему мы должны подходить к величине X1 (т. е. первому прыщу из описанного примера) как к случайной, а не как к постоянной?

Патамушта в теорвере не определено понятие "постоянная величина".
Ergo ваша единица является случайной. Но, чтобы она постоянно присутствовала в каждом испытании, ее - хочешь, не хочешь - надо объявлять достоверным событием. И всю вероятностную массу распределения сбрасывать в точку 1.

Тогда ваша вероятностная модель получаецца самую малость странной: сначала с единичной вероятностью случается (достоверное) событие - некая "постоянная" величина приняла значение, равное 1, а затем к ней еще с некоторой (пуассоновой) вероятностью выпал 0. Так получается пуассонова 1. Для получения пуассоновой двойки надо 1+0+1, etc. Но, увы вам! Вероятность такого "составного" события равна нулю...

Ведь этот первый прыщ есть примененное условие включения пациента в выборку

Ну, если каждый раз отсутствие признака считать пуассоновым нулем... (типа "зуб даю, этот чувак болен прыщами, хотя на морде у него их ровно 0"), то тогда ваша (неслучайная) выборка условна относительно нуля. Это как раз и есть хорошо нам известное zero-truncated Poisson dist. в чистом виде.

Вообще, вся эта комедия объясняется тем, что описанным вами способом никак нельзя модифицировать носитель распределения (0,1,2,...).
Поэтому ребята вынуждены рассматривать ненужные им нули в рамках условного относительно нуля распределения ("нуль-усеченная модель"), трактовать их как число неудач до первого успеха (отрицательное биномиальное распределение и "negative binomial model"), вводить zero-inflated Poisson regression, вспоминать про изученные Кацем в далеком 1963 году Zero-modified Poisson dist. и нет этому конца.

P.S. А вот, a propos, как бы вы перевели словосочетание Zero-inflated Poisson regression? Настоящий пробный камень начинающего переводчика. Или кошмар.

Сообщение отредактировал 100$ - 13.11.2022 - 16:24

[ИНО]

(типа "зуб даю, этот чувак болен прыщами, хотя на морде у него их ровно 0")

Вот-такие не доходят до врача, их еще в дверях заворачивают. Сколько их в популяции, не интересно, да и по количеству прышей у прыщавых едва ли возможно сделать вывод о частоте встречаемости непрыщавых.

Для получения пуассоновой двойки надо 1+0+1, etc. Но, увы вам! Вероятность такого "составного" события равна нулю...

Ту что-то совсем не понял. Двойка у нас никакая не пуассоновская, она получается сложением постоянной 1 с пуассоновской 1. Предположением, что мы не заворачиваем никого в дверях, но за первый "прыщик" на лице пациента считаем... его нос. Он есть у всех и каждого, в отличие от второго и последующего прыщиков, которые действительно прыщики. Врач - инопланетянин, причем очень тупой, и не умеет отличать нос от прыщиков, а может только посчитать общее количество выпуклостей на лице. Вполне себе попадутся ему и с одной, и с двумя и с тремя и т. д., никаких нулевых вероятностей, за исключением собственно нулевых знаний, которых гарантировано не будет (ну, разве то герой Гоголя на прием придет ). По-моему вполне себе адекватная вероятностная модель выходит.

[100$]

Ту что-то совсем не понял. Двойка у нас никакая не пуассоновская, она получается сложением постоянной 1 с пуассоновской 1.

Ну вот, началось в колхозе утро: я ему про вероятности, он мне - про школьную алгебру. Так, чего доброго, и до таблицы умножения дойдет.
А какая же она? Детерминированная, чтоле?

Двойка - самая что ни на есть пуассонова:
- во-первых, несомненно случайная, потому как зависит от случайного второго слагаемого;
- во-вторых, много ли вы знаете дискретных распределений, носителем которых является {0,1,2,3...}?

Если она не пуассоновская, то какая? Из распределения Скеллама, что ли?
А, может быть, Накагами? Трейси-Видома? Уишарта? Уточните, когда не лень.

P.S. Вы распределение-то сгенерируйте: result<-rpois(1000,.685). При параметре распределения Lambda<1 0 - самое вероятное значение. Как от него избавляться? Распределение Пуассона - это же все-таки закон редких событий.

[ИНО]

Что-то я Вас с каждым постом все меньше понимаю. Зачем избавляться от пуассоновского от 0? Он ничуть не мешает, даже наоборот необходим! 1(константа) + 0 (пуассоновский) =1 - как раз самое часто значение в моих подсчетах "прыщиков". И частоты всех последующих значений (1+1=2, 1+2=3...) вполне согласуются с гипотезой о добавлении постоянной единицы к случайной пуассоновской величине. Этот первый прыщик - непременное условие для осмотра больного, достоверное событие, а вот уже количество дополнительных прыщиков - случайная величина. Наезд на алгебру тоже не понял. Прибавление постоянного свободного члена в регрессионной модели Вас тоже смущает?

Сообщение отредактировал ИНО - 14.11.2022 - 10:07

[100$]

Что-то я Вас с каждым постом все меньше понимаю. Зачем избавляться от пуассоновского от 0? Он ничуть не мешает, даже наоборот необходим! 1(константа) + 0 (пуассоновский) =1 - как раз самое часто значение в моих подсчетах "прыщиков". И частоты всех последующих значений (1+1=2, 1+2=3...) вполне согласуются с гипотезой о добавлении постоянной единицы к случайной пуассоновской величине. Этот первый прыщик - непременное условие для осмотра больного, достоверное событие, а вот уже количество дополнительных прыщиков - случайная величина.

Ну, а что тут непонятного? Вы сконструировали некоторое дискретное распределение. Для дискретного распределения можно явно указать вероятность любого события. Вот и укажите вероятность появления суммы вида 1(константа) + 0 (пуассоновский) =1.
Пока этого не сделано, я буду считать это школьной математикой.

Наезд на алгебру тоже не понял. Прибавление постоянного свободного члена в регрессионной модели Вас тоже смущает?

Нет. Потому что в регрессионной модели константа - это такой же (оцениваемый) регрессионный коэффициент, и так же точно является случайной величиной. Для нее так же точно строятся доверительные интервалы, тестируется отличие от нуля.
Так что смущает меня только то, что вы этого не знаете.
Но это не приговор. Люди с этим живут и неплохо себя чувствуют. Вы, например.

Сообщение отредактировал 100$ - 14.11.2022 - 13:04

[ИНО]

Что такое "школьная математика" я не в курсе, равно как и о причине Вашего уничижительного отношения к алгебре.

Вот и укажите вероятность появления суммы вида 1(константа) + 0 (пуассоновский) =1

Та же, что и вероятность 0 в обычной пуассоновской регрессионной модели, построенной по той же выборке, из каждой варианты которой вычтена 1. Формулой это не напишу, увы.

Нет. Потому что в регрессионной модели константа - это такой же (оцениваемый) регрессионный коэффициент, и так же точно является случайной величиной

Может быть оценен, а может быть и взят с неба из теоретических предпосылок, и коэффициенты, кстати, тоже, так что оцениваться будет только ошибка. И это все равно будет статистическая модель.

И да тут все ж форум о прикладной статистике, а не о математической (это разные дисциплины, хотя и связанные). О последней пусть д. ф.-м. н. семиэтажными формулами дискутируют, теоретически обосновывая очередной элегантный чудо-метод имени себя (то, что на практике он может толком не работать, потому что "забыли про овраги" - это уже дело десятое ) . Так что объяснение "на пальцах" для простых смертных типа меня всяко приветствуется.

Сообщение отредактировал ИНО - 15.11.2022 - 09:09

[100$]

Та же, что и вероятность 0 в обычной пуассоновской регрессионной модели, построенной по той же выборке, из каждой варианты которой вычтена 1. Формулой это не напишу, увы.

Ну, если на пальцах, то вы сконструировали некое условное относительно неслучайной величины (единицы) распределение и утверждаете, что в нем условная вероятность равна безусловной.
Хотя это ниоткуда не следует.
Вот, например, если мы захотим сдвинуть стандартное нормальное распределение, мы же не будем городить такую словесную конструкцию "а вот мы дождемся с вероятностью 1 единицу, а потом посмотрим, как реализуется ст.норм.расп." Мы этот сдвиг просто добавим к параметру распределения. И оно будет симметрично теперь уже относительно 1.
Вот бы и с Пуассоном также точно.

[ИНО]

Как обычно, я попробовал развеять теоретический сумрак моделированием. Решил: какая из двух моделей будет иметь более точный прогноз на смоделированных данных, той буду пользоваться на практике. Все ж я биолог, а не математик. Но возник конфуз: прогнозы у сабжа и "моего безграмотного ноу-хау" и их ошибки оказались практически неотличимыми. Более того, при лямбда<1 (а это больше похоже на реальные данные) моя модель работает даже чуточку лучше! Хотя возможно тут дело не непосредственно в лямбде, а в том, что мой алгоритм моделирования усекает объемы выборок, сгенерированных при малых лямбда, но тогда это означает, что расово правильная модель хуже моей ведет себя на малых выборках, так что хрен редьки не слаще.

Уменьшение n на описанную картину не влияет (правда, при n меньше 100 и очень малых лямбда, алгоритм сталкивается с невозможностью оценки модели из-за полного отсутствия ненулевых значений и стопорится, так что в этом диапазоне условий исследование провести не удалось. Наверное, по-хорошему надо применять вместо цикла сопли-функции или как их там, а я в этом не силен. Но для понимания проблемы и имеющееся сойдет.

Итак, прошу указать на мою ошибку, не в теории, которую как тут выяснилось, я постичь не в силах, а именно в моделировании!

Код

# генерация синтетических данных:
lambda<-seq(0.01, 2, 0.01)
simmat<-matrix(NA, 200, 1000)
for(i in 1:200){
lambda_i<-lambda[i]
simmat[i,]<-rpois(1000, lambda_i)}

# предложенная мной раскритикованная выше модель:
prediction<-rep(NA, 200)
MAE<-rep(NA, 200)
for(i in 1:200){
x<-simmat[i,]
df<-as.data.frame(x)
df<-subset(df, x>0)
df$modx<-df$x-1
mod<-glm(modx~1, data=df, family="poisson")
pred<-predict(mod, type="response")+1
prediction[i]<-pred
MAE[i]<-mean(abs(pred-df$x))}

# модель c расово правильно усеченным нулем:
library(VGAM)
prediction2<-rep(NA, 200)
MAE2<-rep(NA, 200)
for(i in 1:200){
x<-simmat[i,]
df<-as.data.frame(x)
df<-subset(df, x>0)
mod<-vglm(x~1, data=df, family="pospoisson")
pred<-predict(mod, type="response")
prediction2[i]<-pred
MAE2[i]<-mean(abs(pred-df$x))}

# пргнозы обеих моделей (моя - красненькая):
plot(lambda, prediction, type="n")
lines(lambda, prediction, col="red")
lines(lambda, prediction2, col="blue")

# средние абсолютные ошибки обеих моделей:
plot(lambda, MAE, type="n")
lines(lambda, MAE, col="red")
lines(lambda, MAE2, col="blue")

# резульаты их сравнения при различных лямбда:
plot(lambda, MAE>MAE2)

Сообщение отредактировал ИНО - 18.11.2022 - 10:49

Эскизы прикрепленных изображений

 

[Диагностик]

при лямбда<1

А что такое лямбда?

[ИНО]

Параметр распределения Пуассона. В моем случае имеется в виду распределение, из которого были сгенерированы синтетические данные. Следует отметить, что варианты с нулевым значением исключены из последующего анализа. По идее, именно исходное распределение, их включающее, пытается реконструировать zero-truncated Poisson regression. Мой же подход работает по-другому: в нем оценивается обычное распределение Пуассона, по выборочным значениям, из которых была вычтена единица. У этого распределения уже иная лямбда, которая здесь не приводится. Да, для простоты в обоих случаях построена регрессия к 1, т, е. просто произведена оценка параметров распределения. Регрессионный анализ как таковой - это уже следующий шаг.

[Диагностик]

Параметр распределения Пуассона.

То есть матожидание. А как здесь ТЕХ использовать?

Сообщение отредактировал Диагностик - 18.11.2022 - 14:32

[100$]

Наконец-то наш естествоиспытатель сам признался, что переливать из пустого в порожнее ему сподручнее безо всякой теории )

Однако теория гласит, что ежели моделировать Пуассоном счетные данные, то с необходимостью надо принять, что эти данные измерены в абсолютной шкале. А в ней, как это знает не только взрослый, но даже карапуз, допустимы только тождественные преобразования, т.е. преобразования вида f(x)=x, переводящие каждый элемент шкалы в самого себя. Ни тебе сложения, ни умножения, ни возведения в степень, ни извлечение квадратного корня. Как из таких данных можно вычесть 1 непонятно.

Далее. В Пуассоново распределение мутирует биномиальное распределение при малой вероятности успеха и неограниченном возрастании N. А что есть биномиальное распределение?
Вот выпало в результате эксперимента из 100 бросков монеты 49 орлов и 51 решка. Как из этих данных вычесть единицу? Это ж фальсификацыя ).

А модель... Ну, что модель: у нее же нет "защиты от дурака". Если она не столкнулась с вычислительными проблемами, то она чего-то там посчитает. Ну и что?