Как построить график функции правдоподобия

Question

Я хочу построить график функции правдоподобия между -pi и pi.

функция логарифмического правдоподобия

llh <- function (teta,x) {

  sum(log((1-cos(x-teta))/(2*pi)))
}

x=c(3.91,4.85,2.28,4.06,3.70,4.04,5.46,3.53,2.28,1.96,2.53,3.88,2.22,3.47,4.82,2.46,2.99,2.54,0.52,2.50)

teta=seq(-4,4, by=0.01)

y = llh(teta,x)

plot(teta, llh(teta,x), pch=16)

Не удалось построить функцию. Вот сообщение об ошибке:

Warning message:
In x - teta :
  longer object length is not a multiple of shorter object length
>   
> plot(teta, llh(teta,x), pch=16)
Error in xy.coords(x, y, xlabel, ylabel, log) : 
  'x' and 'y' lengths differ
In addition: Warning message:
In x - teta :
  longer object length is not a multiple of shorter object length

Answer 1

Как написано, ваша функция будет работать для одного значения teta и несколько x значения или несколько значений teta и один x ценности. В противном случае вы получите неверное значение или предупреждение.

Пример: llh за teta=1 а также teta=2:

> llh(1,x)
[1] -34.88704>
> llh(2,x)
[1] -60.00497

это не то же самое, что:

> llh(c(1,2),x)
[1] -49.50943

И если вы попытаетесь сделать три:

> llh(c(1,2,3),x)
[1] -49.52109
Warning message:
In x - teta :
  longer object length is not a multiple of shorter object length

который в основном происходит от:

> cos(x-c(1,2,3))
[...]
Warning message:
In x - c(1, 2, 3) :
  longer object length is not a multiple of shorter object length

потому что R пытается вычесть вектор длины 3 из вектора длины 20. Он повторяет вектор длины-3 над вектором длины-20 до тех пор, пока это не будет сделано шесть раз, тогда у него останется только два элемента. Хмммм... думает R, я сделаю это, но это выглядит неправильно, поэтому вот предупреждение.

Вы можете либо переписать функцию для работы с векторами для обоих аргументов, либо векторизовать функцию, обернув ее.

> vllh = Vectorize(llh,"teta")
> vllh(c(1,2,3),x)
[1] -34.88704 -60.00497 -67.30765
> plot(teta, vllh(teta,x))

Answer 2

Вам нужно использовать функцию sapply (читай?sapply), так как ваш код не векторизован

plot(teta, sapply(X=teta, FUN=function(teta) llh(teta, x=x)), type="l")

ИЛИ ЖЕ

векторизовать вашу функцию как таковую:

llh2 <- function (teta,x) {
  sapply(X=teta, FUN=function(teta) sum(log((1-cos(x-teta))/(2*pi))) )
}

plot(teta, llh2(teta,x), type="l")

Answer 3

Другое решение с использованием внешнего:

llh <- function (teta,x) {
  X <- outer(x, teta, "-")
  Y <- log((1-cos(X))/(2*pi))
  apply(Y, 2, sum)
}

x=c(3.91,4.85,2.28,4.06,3.70,4.04,5.46,3.53,2.28,1.96,2.53,3.88,2.22,3.47,4.82,2.46,2.99,2.54,0.52,2.50)

teta=seq(-4,4, by=0.01)  

plot(teta, llh(teta,x), pch=16)

Я согласен с разнесенным человеком, что Vectorize интересен!