Rstan на Rstudio MCMC имеет слишком высокое время работы (ограниченное использование доступного ЦП и ОЗУ)

Я новичок в мире Рстана, но мне он очень нужен для моей диссертации. На самом деле я использую сценарий и аналогичный набор данных от парня из Нью-Йоркского университета, который сообщает, что примерное время для аналогичного 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"))

person Mario Luca    schedule 25.04.2015    source источник
comment
Привет, Марио Лука, я пытаюсь запустить тот же код, что и вы (код Пабло Барбера), и у меня такие же трудности. Удалось ли вам это решить?   -  person celacanto    schedule 18.02.2019


Ответы (1)


Во-первых, мои извинения: я бы представил это как комментарий, но у меня недостаточно репутации.

Вот вопрос, который вы задали: «Чем я могу манипулировать, чтобы получить разумный результат за более короткое время?»

Ответ: это зависит от обстоятельств. Вместо того, чтобы представлять вещи в виде двоичной матрицы, пробовали ли вы уменьшить размер матрицы с помощью счетчиков? Основываясь на типе модели, которую вы пытаетесь запустить, я полагаю, что в апостериорной части есть некоторая неидентифицируемость. Не могли бы вы попробовать изменить параметры?

Кроме того, вы можете захотеть запустить CmdStan, если R вызывает проблемы с управлением памятью.

person syclik    schedule 05.09.2015