Предположим, у вас есть случайно сгенерированные переменные в векторе x. Вы можете сделать следующее:

x <- rgeom(1000,0.2)

x_tbl <- table(x)
x_val <- as.numeric(names(x_tbl))
x_df <- data.frame(count=as.numeric(x_tbl), value=x_val)

# Expand to fill in "gaps" in the values caused by 0 counts
all_x_val <- data.frame(value = 0:max(x_val))
x_df <- merge(all_x_val, x_df, by="value", all.x=TRUE)
x_df$count[is.na(x_df$count)] <- 0

# Get theoretical probabilities 
x_df$eprob <- dgeom(x_df$val, 0.2)

# Chi-square test: once with asymptotic dist'n, 
# once with bootstrap evaluation of chi-sq test statistic
chisq.test(x=x_df$count, p=x_df$eprob, rescale.p=TRUE)
chisq.test(x=x_df$count, p=x_df$eprob, rescale.p=TRUE, 
   simulate.p.value=TRUE, B=10000)

[Я думаю, что это было перемещено немного поспешно, несмотря на ваш ответ на вопрос Уубера, поскольку я думаю, что, прежде чем решать проблему «как мне написать этот алгоритм в R», вероятно, более важно разобраться с тем, «что вы выполнение - не лучший подход к проблеме вашей проблемы (которая, безусловно, относится к тому месту, где вы ее разместили). Поскольку он здесь, я буду иметь дело с аспектом «сделать это в R», но я бы настоятельно рекомендовал вам вернуться к вопросу о втором вопросе (как новый пост).]

Во-первых, критерий хи-квадрат немного отличается в зависимости от того, проверяете ли вы

H0: данные получены из геометрического распределения с параметром p

or

H0: данные получены из геометрического распределения с параметром 0,3.

Если вы хотите второй, это довольно просто. Во-первых, с геометрическим, если вы хотите использовать приближение хи-квадрат для распределения тестовой статистики, вам нужно будет сгруппировать соседние ячейки в хвосте. «Обычное» правило - слишком консервативное - предполагает, что вам нужно ожидаемое количество в каждой ячейке не менее 5.

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

Я продолжу, как если бы n было, скажем, 1000 (хотя, если вы тестируете свою генерацию геометрических случайных чисел, это довольно мало).

Во-первых, вычислите ожидаемое количество:

 dgeom(0:20,.3)*1000
 [1] 300.0000000 210.0000000 147.0000000 102.9000000  72.0300000  50.4210000
 [7]  35.2947000  24.7062900  17.2944030  12.1060821   8.4742575   5.9319802
[13]   4.1523862   2.9066703   2.0346692   1.4242685   0.9969879   0.6978915
[19]   0.4885241   0.3419669   0.2393768

Предупреждение, dgeom и друзья идут от x = 0, а не x = 1; хотя вы можете перенести входы и выходы на функции R, это будет намного проще, если вы вычтете 1 из всех своих геометрических значений и проверите это. Я буду действовать так, как если бы из вашей выборки вычли 1, так что она идет от 0.

Я отключу это на 15-м семестре (x = 14) и сгруппирую 15+ в отдельную группу (в данном случае - единственную группу). Если вы хотите следовать эмпирическому правилу «больше пяти», вы бы отключили его после 12-го члена (x = 11). В некоторых случаях (например, при меньшем p) вам может потребоваться разделить хвост на несколько ячеек, а не на один.

> expec <- dgeom(0:14,.3)*1000
> expec <- c(expec, 1000-sum(expec))
> expec
 [1] 300.000000 210.000000 147.000000 102.900000  72.030000  50.421000
 [7]  35.294700  24.706290  17.294403  12.106082   8.474257   5.931980
[13]   4.152386   2.906670   2.034669   4.747562

Последняя ячейка - это категория «15+». Нам также нужны вероятности.

Сейчас у нас еще нет образца; Я просто сгенерирую один:

y <- rgeom(1000,0.3)

но теперь нам нужна таблица наблюдаемых подсчетов:

 (x <- table(factor(y,levels=0:14),exclude=NULL))

   0    1    2    3    4    5    6    7    8    9   10   11   12   13   14 <NA> 
 292  203  150   96   79   59   47   25   16   10    6    7    0    2    5    3 

Теперь вы можете напрямую вычислить хи-квадрат, а затем вычислить значение p:

> (chisqstat <- sum((x-expec)^2/expec))
[1] 17.76835
(pval <- pchisq(chisqstat,15,lower.tail=FALSE))
[1] 0.2750401

но вы также можете заставить R сделать это:

> chisq.test(x,p=expec/1000)

        Chi-squared test for given probabilities

data:  x 
X-squared = 17.7683, df = 15, p-value = 0.275

Warning message:
In chisq.test(x, p = expec/1000) :
  Chi-squared approximation may be incorrect

Теперь случай для неопределенного p аналогичен, но (насколько мне известно) вы больше не можете заставить chisq.test делать это напрямую, вам нужно сделать это первым способом, но вы должны оценить параметр на основе данных (по максимальной вероятности или минимальный хи-квадрат), а затем проверьте, как указано выше, но у вас на одну степень свободы меньше для оценки параметра.

См. Пример выполнения хи-квадрат для Пуассона с предполагаемым параметром здесь; геометрия следует тому же подходу, что и выше, с настройками, как в ссылке (работа с неизвестным параметром, включая потерю 1 степени свободы).

В пакете vcd есть функция goodfit, описанная как Goodness-of-fit Tests для дискретных данных.

G.fit <- goodfit(x, type = "nbinomial", par = list(size = 1))

Я собирался использовать код, который вы разместили в предыдущем вопросе, но теперь оказалось, что вы удалили этот код. Я считаю это оскорблением. Вы используете этот форум, чтобы собрать ответы на домашнее задание, а затем стереть его, чтобы удалить улики? (Удаленные вопросы все еще могут быть видны тем из нас, у кого достаточно репутации, а интерфейс предотвращает удаление вопроса с одобренными ответами, поэтому вы не сможете удалить этот.)

Создайте график QQ для тестирования геометрически распределенного образца

sim.geometric <- function(nvals)
{
    p <- 0.3
    u <- runif(nvals)
    ceiling(log(u)/log(1-p))
}

png("QQ.png")
qqplot(qgeom(ppoints(100),prob=0.3), sim.res,
       main = expression("Q-Q plot for" ~~ {G}[n == 100]))
dev.off()

sim.res <- sim.geometric(500)
qqplot(jitter(qgeom(ppoints(500),prob=0.3)), jitter(sim.res),
       main = expression("Q-Q plot for" ~~ {G}[n == 100]), ylim=c(0,max( qgeom(ppoints(500),prob=0.3),sim.res )),
xlim=c(0,max( qgeom(ppoints(500),prob=0.3),sim.res )))
 qqline(sim.res, distribution = function(p) qgeom(p, 0.3),
       prob = c(0.25, 0.75), col = "red")