байесовская PCA с использованием PyMC

Я пытаюсь реализовать байесовский PCA, используя библиотеку PyMC для python. Но я застрял там, где определяю координаты более низкого измерения...

Модель

x = Wz + e

где x — вектор наблюдения, W — матрица преобразования, а z — вектор координат меньшего размера.

Сначала я определяю распределение для матрицы преобразования W. Каждый столбец взят из нормального распределения (для простоты нулевое среднее значение и тождественная ковариация).

def W_logp(value):
   logLikes = np.array([multivariate_normal.logpdf(value[:,i], mean=np.zeros(dimX), cov=1) for i in range(0, dimZ)])
   return logLikes.sum()

def W_random():
   W = np.zeros([dimX, dimZ])
   for i in range(0, dimZ):
      W[:,i] = multivariate_normal.rvs(mean=np.zeros(dimX), cov=1)
   return W

w0 = np.random.randn(dimX, dimZ)

W = pymc.Stochastic(
   logp = W_logp,
   doc = 'Transformation',
   name = 'W',
   parents = {},
   random = W_random,
   trace = True,
   value = w0,
   dtype = float,
   rseed = 116.,
   observed = False,
   cache_depth = 2,
   plot = False,
   verbose = 0)

Затем я хочу определить распределение для z, которое снова является многомерным нормальным (нулевое среднее и тождественная ковариация). Однако мне нужно нарисовать z для каждого наблюдения отдельно, а W общее для всех. Итак, я попытался

z = pymc.MvNormal('z', np.zeros(dimZ), np.eye(dimZ), size=N)

Однако pymc.MvNormal не имеет параметра size. Так что выдает ошибку. Следующим шагом будет

m = Data.mean(axis=0) + np.dot(W, z)
obs = pymc.MvNormal('Obs', m, C, value=Data, observed=True)

Я не приводил спецификацию для C выше, так как сейчас это не имеет значения. Есть идеи как реализовать?

Спасибо

РЕДАКТИРОВАТЬ

После ответа Криса Фоннесбека я изменил свой код следующим образом

numD, dimX = Data.shape
dimZ = 3
mm = Data.mean(axis=0)

tau = pymc.Gamma('tau', alpha=10, beta=2)
tauW = pymc.Gamma('tauW', alpha=20, beta=2, size=dimZ)

@pymc.deterministic(dtype=float)
def C(tau=tau):
    return (tau)*np.eye(dimX)

@pymc.deterministic(dtype=float)
def CW(tau=tauW):
    return np.diag(tau)

W = [pymc.MvNormal('W%i'%i, np.zeros(dimZ), CW) for i in range(dimX)]
z = [pymc.MvNormal('z%i'%i, np.zeros(dimZ), np.eye(dimZ)) for i in range(numD)]
mu = [pymc.Lambda('mu%i'%i, lambda W=W, z=z: mm + np.dot(np.array(W), np.array(z[i]))) for i in range(numD)]

obs = [pymc.MvNormal('Obs%i'%i, mu[i], C, value=Data[i,:], observed=True) for i in range(numD)]

model = pymc.Model([tau, tauW] + obs + W + z)
mcmc = pymc.MCMC(model)

Но на этот раз он пытается выделить огромное количество памяти (более 8 ГБ) при запуске pymc.MCMC(model) с numD=45 и dimX=504. Даже когда я пробую это только с numD=1 (то есть создаю только 1 z, mu и obs), он делает то же самое. Есть идеи, почему?


person ahmethungari    schedule 23.10.2014    source источник


Ответы (2)


К сожалению, PyMC не позволяет легко определять векторы многомерного стохастика. Надеюсь, мы сможем сделать это в PyMC 3. Сейчас вам нужно указать это с помощью контейнера. Например:

z = [pymc.MvNormal('z_%i' % i, np.zeros(dimZ), np.eye(dimZ)) for i in range(N)]
person Chris Fonnesbeck    schedule 28.10.2014
comment
Спасибо за ответ. Позвольте мне добавить еще один, есть ли способ использовать таблички как в графических моделях, в PyMC? Действительно, было бы очень полезно иметь какой-либо пример, сначала дающий графическую модель, а затем кодирующий в PyMC. - person ahmethungari; 29.10.2014
comment
Пожалуйста, смотрите мой РЕДАКТИРОВАТЬ в вопросе - person ahmethungari; 29.10.2014

Что касается проблемы с памятью, попробуйте использовать другой сервер для трассировки. По умолчанию ("ram") все хранится в оперативной памяти. Вместо этого вы можете попробовать что-то вроде "pickle" или "sqlite".

Что касается обозначения номерного знака, возможно, это то, что мы могли бы реализовать для PyMC 3. Не стесняйтесь создавать проблему с предложением об этом в наш трекер ошибок.

person Chris Fonnesbeck    schedule 30.10.2014