Отбор проб для генерации нормальных образцов из образцов Коши - PullRequest
0 голосов
/ 06 июля 2018

Я попытал счастья при кодировании метода выборки с отбраковкой, чтобы создать выборку, которая соответствует нормальному распределению. Образцы выглядят как нормальные распределения на первый взгляд, но значение p теста Шапиро-Уилка всегда <0,05. Я действительно не знаю, где я ошибся, и я получил только псевдокод от своего учителя (это НЕ домашняя работа). Любая помощь приветствуется. Ниже мой код: </p>

f <- function(x,m,v) {    #target distribution, m=mean,v=variance
  dnorm(x,m,sqrt(v))
}

g <- function(x,x0,lambda) {  #cauchy distribution for sampling
  dcauchy(x,x0,lambda)
}

genSamp <- function(n,m,v) {  #I want the user to be able to choose mean and sd
                              #and size of the sample
  stProbe <- rep(0,n)         #the sample vector
  interval = c(m-10*sqrt(v),m+10*sqrt(v)) #wanted to go sure that everything
                                          #is covered, so I took a range
                                          #that depends on the mean
  M = max(f(interval,m,v)/g(interval,m,v))  #rescaling coefficient, so the cauchy distribution
                              #is never under the normal distribution
  #I chose x0 = m and lambda = v, so the cauchy distribution is close to a
  #the target normal distribution

  for (i in 1:n) {
    repeat{
      x <- rcauchy(1,m,v)
      u <- runif(1,0,max(f(interval,m,v)))
      if(u < (f(x,m,v)/(M*g(x,m,v)))) {
        break
      }
    }
    stProbe[i] <- x
  }

  return(stProbe)
}

Тогда я попробовал это с:

test <- genSamp(100,2,0.5)
hist(test,prob=T,breaks=30)#looked not bad
shapiro.test(test) #p-value way below 0.05

Заранее благодарю за помощь.

Ответы [ 2 ]

0 голосов
/ 06 июля 2018

На самом деле, первое, что я проверил, это среднее значение выборки и дисперсия выборки. Когда я рисую 1000 выборок с вашим genSamp, я получаю среднее значение выборки на уровне 2, но дисперсию выборки на уровне 2,64, что далеко от целевого значения 0,5.

Первая проблема связана с вашим вычислением M. Обратите внимание, что:

interval = c(m - 10 * sqrt(v), m + 10 * sqrt(v))

дает только 2 значения, а не сетку из одинаково расположенных точек на интервале. На 10 стандартных отклонений от среднего значения нормальная плотность равна почти 0, поэтому M составляет почти 0. Вам нужно сделать что-то вроде

interval <- seq(m - 10 * sqrt(v), m + 10 * sqrt(v), by = 0.01)

Вторая проблема - генерация единой случайной величины в вашем repeat. Почему вы делаете

u <- runif(1,0,max(f(interval,m,v)))

Вы хотите

u <- runif(1, 0, 1)

С этими исправлениями я проверил, что genSamp получает правильное среднее значение выборки и дисперсию выборки. Образцы проходят как тест Шапиро – Вилка, так и тест Колмогорова-Смирнова (?ks.test).


Полный рабочий код

f <- function(x,m,v) dnorm(x,m,sqrt(v))

g <- function(x,x0,lambda) dcauchy(x,x0,lambda)

genSamp <- function(n,m,v) {

  stProbe <- rep(0,n)
  interval <- seq(m - 10 * sqrt(v), m + 10 * sqrt(v), by = 0.01)
  M = max(f(interval,m,v)/g(interval,m,v))

  for (i in 1:n) {
    repeat{
      x <- rcauchy(1,m,v)
      u <- runif(1,0,1)
      if(u < (f(x,m,v)/(M*g(x,m,v)))) break
      }
    stProbe[i] <- x
    }

  return(stProbe)
  }

set.seed(0)
test <- genSamp(1000, 2, 0.5)
shapiro.test(test)$p.value
#[1] 0.1563038

ks.test(test, rnorm(1000, 2, sqrt(0.5)))$p.value
#[1] 0.7590978
0 голосов
/ 06 июля 2018

У вас есть

f <- function(x,m,v) {    #target distribution, m=mean,v=variance
  dnorm(x,e,sqrt(v))
}

какие выборки имеют среднее значение e, но это никогда не определяется.

...