Gostaria de desenhar com eficiência amostras de sujeitas à restrição de que .
Gostaria de desenhar com eficiência amostras de sujeitas à restrição de que .
A resolução formal desse problema primeiro requer uma definição adequada de um
" distribuição sujeita à restrição de que "
A maneira natural é definir a distribuição de condicional a . E para aplicar isso condicional ao caso . Se usarmos coordenadas polares , o jacobiano da transformação é Portanto, a densidade condicional da distribuição de| | X | | = ϱ ϱ = 1 x 1 = ϱ cos ( θ 1 ) ϱd-1
Conclusão: Essa densidade difere da simples aplicação da densidade Normal a um ponto na esfera unitária por causa do jacobiano.
O segundo passo é considerar a densidade de destino e projete um algoritmo de Monte Carlo da cadeia de Markov para explorar o espaço de parâmetros . Minha primeira tentativa seria em um amostrador de Gibbs, inicializado no ponto na esfera mais próxima de , ou seja,e procedendo um ângulo de cada vez, da maneira Metropolis-in-Gibbs: [0,π] d - 2 ×[0,2π]μμ
As escalas , , , podem ser dimensionadas com base nas taxas de aceitação das etapas, na direção de uma meta ideal de .
Aqui está um código R para ilustrar o acima, com valores padrão para e :Σ
library(mvtnorm)
d=4
target=function(the,mu=1:d,sigma=diag(1/(1:d))){
carte=cos(the[1])
for (i in 2:(d-1))
carte=c(carte,prod(sin(the[1:(i-1)]))*cos(the[i]))
carte=c(carte,prod(sin(the[1:(d-1)])))
prod(sin(the)^((d-2):0))*dmvnorm(carte,mean=mu,sigma=sigma)}
#Gibbs
T=1e4
#starting point
mu=(1:d)
mup=mu/sqrt(sum(mu^2))
mut=acos(mup[1])
for (i in 2:(d-1))
mut=c(mut,acos(mup[i]/prod(sin(mut))))
thes=matrix(mut,nrow=T,ncol=d-1,byrow=TRUE)
delta=rep(pi/2,d-1) #scale
past=target(thes[1,]) #current target
for (t in 2:T){
thes[t,]=thes[t-1,]
for (j in 1:(d-1)){
prop=thes[t,]
prop[j]=prop[j]+runif(1,-delta[j],delta[j])
prop[j]=prop[j]%%(2*pi-(j<d-1)*pi)
prof=target(prop)
if (runif(1)<prof/past){
past=prof;thes[t,]=prop}
}
}
x E [ ( x - μ ) 2 ] ˜ = 1 não é estritamente possível, pois é uma variável aleatória (contínua). Se você deseja que ele varie 1, ou seja, (onde o til significa que estimamos a variação), você precisará exigir que sua variação seja . No entanto, essa demanda pode entrar em conflito com . Ou seja, para obter amostras com essa variação, você precisa que a diagonal de seja igual a . 1 ΣΣ1
Para fazer o exemplo dessa distribuição em geral, você pode gerar normais padrão iid e multiplicar por , a raiz quadrada de e adicionar os meios . Σ u