m**c
发帖数: 88 | 1 有个问题一直没解决,找不到人问,发上来希望能有高手指点一下.很简单的模型,线性高
斯模型, Y = A0 + A*X + e, 已知观测向量Y(5维的),系数矩阵 (A0,A)以及方差矩阵(V
),估计X(7维的). 现在的问题是,X维数太高,采样不收敛. 降低X的维数,比如我只估计X
[1],X[2],其他X分量的先验分布的方差设置的很小,这样相当于其他分量为固定值了,这
时后采样效果较好,算法收敛,但是进一步提高X的维数,比如下面这段程序里,估计X[1],
X[2],X[4],X[5],X[7],结果很不好.
想问问大家我这个程序这样写有没有问题?有没有办法改进,使得估计高维的X是结果仍
然可以接受?
model{
for (i in 1:n){
y[i,1:5]~ dmnorm( mu[], T[,] )
}
# computing the mu
for (i in 1:P){
mu[i] <- inprod(A[i,],x[,])+ A0[i,1]
}
... 阅读全帖 |
|
发帖数: 1 | 2 here's my code:
bayes.mod<-function(){
#proficiency model
for (i in 1:10){
for (k in 1:4){
pi[i,k]<-exp(lambda0[i]+lambda1[i]*theta[i])/(1+exp(lambda0[i]+lambda1
[i]*theta[i]))
alpha[i,k]~dcat(pi[i,k]) #categorical density
}
lambda0[i]~dunif(-1,1)
lambda1[i]~dunif(0,2)
theta[i]~dnorm(0,1)}
#evidence model
for ( j in 1:10){
for (k in 1:4){
gamma[j,k] <- rr[j,k]*q[j,k]
rr [j,k]~ dnorm(1,3)
q<-q
}
beta[j]~dnorm(mu,1)
}
mu~d... 阅读全帖 |
|
W*****r
发帖数: 193 | 3 model {
for (j in 1:J){
y[j] ~ dnorm(theta[j], tau.y[j])
theta[j] ~dnorm(mu.theta, tau.theta)
tau.y[j] <- pow(sigma.y[j], -2)
}
mu.theta ~ dnorm(0, 1.0E-6)
tau.theta <- pow(sigma.theta, -2)
sigma.theta ~ dunif (0, 1000)
either<-max(theta1,theta2，theta3,theta4, theta5, theta6, the7, theta8)
问题在这一行，我共有8个theta值。这样输入说“没有定义theta1...“
我想请教一下：应该怎么输入？
多谢多谢。会发包子的。
}
DATA
list(J=8, y = c(28, 8, -3, 7, -1,1,18,12), sigma.y=c(15,10,16,11,9,11,
10,18))
INITIAL VALUES
list(theta= c(0,0,0,0,0,0,0,0), mu.th... 阅读全帖 |
|
W*****r
发帖数: 193 | 4 model {
for (j in 1:J){
y[j] ~ dnorm(theta[j], tau.y[j])
theta[j] ~dnorm(mu.theta, tau.theta)
tau.y[j] <- pow(sigma.y[j], -2)
}
mu.theta ~ dnorm(0, 1.0E-6)
tau.theta <- pow(sigma.theta, -2)
sigma.theta ~ dunif (0, 1000)
either<-max(sort[theta[1], theta[2],theta[3],theta[4]],sort[theta[5],
theta[6],theta[7],theta[8]]))
}
我又试了一下，说有variable未定义。BTW，我只用两个theta值的时候，either可以
set成node，density也可以显示。
多谢多谢。
DATA
list(J=8, y = c(28, 8, -3, 7, -1,1,18,12), sigma.y=c(15,10,16,11,9,11,
10,18))
INITIAL VALUES... 阅读全帖 |
|
m**c
发帖数: 88 | 5 model{
for (i in 1:n){
y[i,1:2]~ dmnorm( mu[], T[,] )
}
for (i in 1:P){
mu[i] <- inprod(A[i,],x[,])+ A0[i,1]
}
x[1,1] ~ dnorm(0.0, 1.0E-4)
x[2,1] ~ dnorm(0.0, 1.0E-4)
#x[1,1] ~ dunif(3,7)
#x[2,1] ~ dunif(4,8)
T[1:P,1:P] <- inverse(V[,])
}
INITS
list(x = structure(.Data = c(3,6), .Dim =c(2,1)))
DATA (RECT.)
list( n=10, P=2,
A = structure(.Data = c(1.5,1.2, 2, 0.7),.Dim =c(2,2)),
A0 = structure(.Data... 阅读全帖 |
|
z**********9
发帖数: 133 | 6 为什么出现这个error：
Error in jags.model(textConnection(model_string), data = list(Y = Y, theta =
theta, :
Error parsing model file:
syntax error on line 5 near "("
mu0 <- 0
s20 <- 1000
a <- 0.01
b <- 0.01
#
library(rjags)
model_string <- "model{
# Likelihood
for(i in 1:n){
Y[i] ~ dpois(lamda(theta[i]))
lamda(theta[i]) <- alpha_0 + sum(beta[]*exp(-0.5(theta[i]-gamma[])^2
/tau[]^2))
}
# Defining Priors
alpha_0 ~ dnorm(0, 0.0001)
for(j in 1:6){
beta... 阅读全帖 |
|
s**u
发帖数: 383 | 7 刚开始学WINBUG, 写了个小程序, model 问题, 但是load data 时说 expected
variable name. 请问哪里出错了. 谢谢.
model {
for (i in 1:9) {n[i] <- nwin1[i] + nwin2[i]
winn1[i] ~dbin(p[i],n[i])
logit(p[i]) <- d[dind1[i]]-d[dind2[i]]}
for ( i in 1:2) {d[i] ~ dnorm(0, 0.001)}
for ( i in 4:5) {d[i] ~ dnorm(0, 0.001)}
d[3] <-0
}
list (dind1 = c(1,1,1,1,2,2,2,3,4),
dind2 = c(2,3,4,5,3,4,5,5,5),
nwin1 = c(9,12,6,27,9,12,12,2,2),
nwin2 = c(10,2,3,2,5,3,3,0,4))
list(d =c(0, 0, NA,0,0)) |
|
s*******1
发帖数: 146 | 8 P(X>3)=dnorm(3, mean=1, sd=2)=0.121;
P(X>2)=dnorm(2, mean=1, sd=2)=0.176
P(X>2)-P(X>3)=0.055 |
|
l****a
发帖数: 336 | 9 模型通过, data load不上去, 总是提示"Expected variable name". 我查了code， 没
看到未定义参数阿。
跪求大侠指点问题在哪儿？
model{
#specify model
for(i in 1:n){
y[i] ~ dnorm(mu[i], 1) #likelihood
mu[i] <- inprod(wPrime[i,], eta[ ]) #process model
}
#prior distribution of coefficients
for(k in 1:R){
eta[k] ~ dnorm(0, 1)}
} # end of model
#data
list(
R=3, n=5,
y=c(0.8295724, 1.0092876, 0.6525292, 0.6990704, 0.3756401),
wPrime = structure(
.Data = c(
... 阅读全帖 |
|
m**c
发帖数: 88 | 10 model{
for (t in 1:n){
y[t,1:5]~ dmnorm( mu[], T[,] )
}
index ~ dcat(Prob[])
for (i in 1:H) {
x_spike[i] ~ dnorm(1,1.0E10)
x_slap[i] ~ dnorm(0.98,100)
lambda[i] <- Index[index,i]
x[i,1] <- (1- lambda[i])* x_spike[i] + lambda[i]* x_slap[i]
}
for (j in 1:P){
mu[j] <- inprod(A[j,],x[,])+ A0[j,1]
}
T[1:P,1:P] <- inverse(V[,])
}
DAtA
.................
.................
.................
END
这个WinBUGS程序我可以在Sampl... 阅读全帖 |
|
|
f***a
发帖数: 329 | 12 luo tie chu lai ba ...
##### Cell growth program for Nested Sampling ####
#In 2-D parameter space, expand from a starting point
#until all points on boundary reach threshold likelihood.
###################
rm(list=ls(all=T))
###################
###### Functions (run once)#######
LL <- function(mu,sigma){sum(log(dnorm(dat,mu,sigma)))}
#log-likelihood function
walk <- function(p){
p.x <- p.nx <- p.y <- p.ny <- NULL
if(p$x==0) p.x <- c(x=0,y=0,nx=1,ny=0,a=p$a+da,b=p$b)
if(p$nx==0) p.nx <- c(x=1,y=0 |
|
a*********e
发帖数: 1233 | 13 本人初用WinBUGS。
假设有模型 mu[i] <- beta0+ beta1*p[i]
其中i=100，
mu[i] 是100个已知数据，
p[i]是生成的数据。
问题是，p[i]分成两部分，例如前50个是均匀分布，后50个是正态分布。
下面的代码显然不能运行
for (i in 1:m){
y[i] = dunif(0,10)
}
for(i in m+1:n){
y[i] = dnorm(0,1)
}
m=50, n=100
请问有什么方法能生成这样的数据，或者说能否把生成的几个数据合并成一个变量。
谢谢。 |
|
n******0
发帖数: 61 | 14 google some introduction or read R extesion,
example cited from "An Introduction to the .C Interface to R":
#include
#include
void kernel_smooth(double *x, int *n, double *xpts, int *nxpts,
double *h, double *result)
{
int i, j;
double d, ksum;
for(i=0; i < *nxpts; i++) {
ksum = 0;
for(j=0; j < *n; j++) {
d = xpts[i] - x[j];
ksum += dnorm(d / *h, 0, 1, 0); // see use R function in C code
}
result[i] = ksum / ((*n) * (*h));
}
} |
|
A**P
发帖数: 260 | 15 在WINBUGS下面写了这么简单的一段，model syntax通过了，可是data就是不能load。
试过list和matrix,都不行。请高手看看为什么。多谢。
model {
#model
nI ~ dbin(p, N)
nS ~ dbin(s,nI)
lor <- log(p*s/(1-p*s))
orm ~ dnorm(lor, orp)
#priors
p ~ dbeta(1,1)
s ~ dbeta(1,1)
}
#Data
list (nI=500, nS=200, N=1421, orm= -1.7, orp= 100)
nI [] nS [] N [] orm [] orp []
500 200 1421 -1.7 100
END |
|
