Я новичок в мире Рстана, но мне он очень нужен для моей диссертации. На самом деле я использую сценарий и аналогичный набор данных от парня из Нью-Йоркского университета, который сообщает, что примерное время для аналогичного DS составляет около 18 часов. Однако, когда я пытаюсь запустить свою модель, она не сделает больше 10% за 18 часов. Таким образом, я прошу небольшой помощи, чтобы понять, что я делаю не так и как повысить эффективность.
Я использую модель двух цепочек на 500 и 100 разминок с функцией Bernoulli_logit по 5 параметрам, пытаясь оценить 2 из них с помощью процедуры MC без разворота. (на каждом шаге он извлекает из случайной нормы каждый параметр, затем оценивает y и сравнивает его с фактическими данными, чтобы увидеть, подходят ли новые параметры к данным лучше)
y[n] ~ bernoulli_logit( alpha[kk[n]] + beta[jj[n]] - gamma * square( theta[jj[n]] - phi[kk[n]] ) );
(n составляет около 10 млн.) Мои данные представляют собой матрицу размером 10.000x1004 из нулей и единиц. Подводя итог, это матрица о людях, подписанных на политиков в твиттере, и я хочу оценить их политические идеи с учетом того, за кем они подписаны. Я запускаю модель на RStudio с R x64 3.1.1 на Win8 Professional, 6 бит, четырехъядерный процессор I7 с оперативной памятью 16 ГБ. Проверяя производительность, rsession использует не более 14% ЦП и 6 ГБ оперативной памяти, хотя еще 7 ГБ бесплатны. Пытаясь выполнить субдискретизацию до матрицы размером 10.000x250, я заметил, что вместо этого будет использоваться менее 1,5 ГБ. Однако я попробовал эту процедуру с набором данных 50x50, и она сработала отлично, поэтому в процедуре нет ошибки. Rsession открывает 8 потоков, я вижу активность на каждом ядре, но ни один из них не занят полностью. Интересно, почему мой компьютер не работает в меру своих возможностей, и может ли быть какое-то узкое место, ограничение или настройка, которые мешают ему это сделать. R - 64-битный (только что проверен), поэтому Rstan должен быть (хотя у меня были некоторые трудности при установке, и это могло испортить некоторые параметры)
вот что происходит, когда я его компилирую
Iteration: 1 / 1 [100%] (Sampling)
# Elapsed Time: 0 seconds (Warm-up)
# 11.451 seconds (Sampling)
# 11.451 seconds (Total)
SAMPLING FOR MODEL 'stan.code' NOW (CHAIN 2).
Iteration: 1 / 1 [100%] (Sampling)
# Elapsed Time: 0 seconds (Warm-up)
# 12.354 seconds (Sampling)
# 12.354 seconds (Total)
в то время как, когда я запускаю его, он просто работает часами, но никогда не выходит за пределы 10% первой цепочки (в основном потому, что я прервал его после того, как мой компьютер собирался расплавиться).
Iteration: 1 / 500 [ 0%] (Warmup)
и имеет этот параметр:
stan.model <- stan(model_code=stan.code, data = stan.data, init=inits, iter=1, warmup=0, chains=2)
## running modle
stan.fit <- stan(fit=stan.model, data = stan.data, iter=500, warmup=100, chains=2, thin=thin, init=inits)
пожалуйста, помогите мне найти, что замедляет процедуру (и если ничего не происходит, чем я могу манипулировать, чтобы получить разумный результат за более короткое время?).
Я благодарю вас заранее,
ML
вот модель (от Пабло Барбера, Нью-Йоркский университет)
n.iter <- 500
n.warmup <- 100
thin <- 2 ## this will give up to 200 effective samples for each chain and par
Adjmatrix <- read.csv("D:/TheMatrix/Adjmatrix_1004by10000_20150424.txt", header=FALSE)
##10.000x1004 matrix of {0, 1} with the relationship "user i follows politician j"
StartPhi <- read.csv("D:/TheMatrix/StartPhi_20150424.txt", header=FALSE)
##1004 vector of values [-1, 1] that should be a good prior for the Phi I want to estimate
start.phi<-ba<-c(do.call("cbind",StartPhi))
y<-Adjmatrix
J <- dim(y)[1]
K <- dim(y)[2]
N <- J * K
jj <- rep(1:J, times=K)
kk <- rep(1:K, each=J)
stan.data <- list(J=J, K=K, N=N, jj=jj, kk=kk, y=c(as.matrix(y)))
## rest of starting values
colK <- colSums(y)
rowJ <- rowSums(y)
normalize <- function(x){ (x-mean(x))/sd(x) }
inits <- rep(list(list(alpha=normalize(log(colK+0.0001)),
beta=normalize(log(rowJ+0.0001)),
theta=rnorm(J), phi=start.phi,mu_beta=0, sigma_beta=1,
gamma=abs(rnorm(1)), mu_phi=0, sigma_phi=1, sigma_alpha=1)),2)
##alpha and beta are the popularity of the politician j and the propensity to follow people of user i;
##phi and theta are the position on the political spectrum of pol j and user i; phi has a prior given by expert surveys
##gamma is just a weight on the importance of political closeness
library(rstan)
stan.code <- '
data {
int<lower=1> J; // number of twitter users
int<lower=1> K; // number of elite twitter accounts
int<lower=1> N; // N = J x K
int<lower=1,upper=J> jj[N]; // twitter user for observation n
int<lower=1,upper=K> kk[N]; // elite account for observation n
int<lower=0,upper=1> y[N]; // dummy if user i follows elite j
}
parameters {
vector[K] alpha;
vector[K] phi;
vector[J] theta;
vector[J] beta;
real mu_beta;
real<lower=0.1> sigma_beta;
real mu_phi;
real<lower=0.1> sigma_phi;
real<lower=0.1> sigma_alpha;
real gamma;
}
model {
alpha ~ normal(0, sigma_alpha);
beta ~ normal(mu_beta, sigma_beta);
phi ~ normal(mu_phi, sigma_phi);
theta ~ normal(0, 1);
for (n in 1:N)
y[n] ~ bernoulli_logit( alpha[kk[n]] + beta[jj[n]] -
gamma * square( theta[jj[n]] - phi[kk[n]] ) );
}
'
## compiling model
stan.model <- stan(model_code=stan.code,
data = stan.data, init=inits, iter=1, warmup=0, chains=2)
## running modle
stan.fit <- stan(fit=stan.model, data = stan.data,
iter=n.iter, warmup=n.warmup, chains=2,
thin=thin, init=inits)
samples <- extract(stan.fit, pars=c("alpha", "phi", "gamma", "mu_beta",
"sigma_beta", "sigma_alpha"))