优化 TensorFlow 2 中许多伪数据实现的功能
Optimise function for many pseudodata realisations in TensorFlow 2
我的最终目标是模拟似然比检验统计，然而，我遇到的核心问题是我不明白如何让 TensorFlow 2 对不同的数据输入执行许多优化。这是我的尝试，希望它能让您了解我正在尝试的事情：
import tensorflow as tf
import tensorflow_probability as tfp
from tensorflow_probability import distributions as tfd
import numpy as np
# Bunch of independent Poisson distributions that we want to combine
poises0 = [tfp.distributions.Poisson(rate = 10) for i in range(5)]
# Construct joint distributions
joint0 = tfd.JointDistributionSequential(poises0)
# Generate samples
N = int(1e3)
samples0 = joint0.sample(N)
# Now we need the same distributions but with floating parameters,
# and need to define the function to be minimised
mus = [tf.Variable(np.random.randn(), name='mu{0}'.format(i)) for i in range(5)]
#@tf.function
def loss():
poises_free = [tfp.distributions.Poisson(rate = mus[i]) for i in range(5)]
joint_free = tfd.JointDistributionSequential(poises_free)
# Construct (half of) test statistic
return -2*(joint_free.log_prob(samples0))
# Minimise (for all samples? Apparently not?)
opt = tf.optimizers.SGD(0.1).minimize(loss,var_list=mus)
print(mus)
print(loss())
print(opt)
quit()
输出：
[<tf.Variable 'mu0:0' shape=() dtype=float32, numpy=53387.016>, <tf.Variable 'mu1:0' shape=() dtype=float32, numpy=2540.568>, <tf.Variable 'mu2:0' shape=() dtype=float32, numpy=-5136.6226>, <tf.Variable 'mu3:0' shape=() dtype=float32, numpy=-3714.5227>, <tf.Variable 'mu4:0' shape=() dtype=float32, numpy=1062.9396>]
tf.Tensor(
[nan nan nan nan ... nan nan nan], shape=(1000,), dtype=float32)
<tf.Variable 'UnreadVariable' shape=() dtype=int64, numpy=1>
最后我想计算检验统计量
q = -2*joint0.log_prob(samples0) - loss()
并表明它服从具有 5 个自由度的卡方分布。
我是 TensorFlow 的新手，所以也许我做的完全错了，但我希望你明白我想要什么。
编辑：
所以我又玩了一会儿，我想 TensorFlow 根本不会像我假设的那样对输入张量进行并行优化。或者也许可以，但我需要进行不同的设置，即可能一次给它一个输入参数的张量和一个巨大的联合损失函数来进行所有最小化？
我也试过用一个简单的循环来做事情，看看会发生什么。正如预测的那样，它慢得可怜，但我什至没有得到正确的答案：
poises0 = [tfp.distributions.Poisson(rate = 10) for i in range(5)]
joint0 = tfd.JointDistributionSequential(poises0)
N = int(5e2)
samples0 = joint0.sample(N)
mus = [tf.Variable(10., name='mu{0}'.format(i)) for i in range(5)]
#@tf.function
def loss(xi):
def loss_inner():
poises_free = [tfp.distributions.Poisson(rate = mus[i]) for i in range(5)]
joint_free = tfd.JointDistributionSequential(poises_free)
# Construct (half of) test statistic
return -2*(joint_free.log_prob(xi))
return loss_inner
# Minimise
# I think I have to loop over the samples... bit lame. Can perhaps parallelise though.
q = []
for i in range(N):
xi = [x[i] for x in samples0]
opt = tf.optimizers.SGD(0.1).minimize(loss=loss(xi),var_list=mus)
q += [-2*joint0.log_prob(xi) - loss(xi)()]
fig = plt.figure()
ax = fig.add_subplot(111)
sns.distplot(q, kde=False, ax=ax, norm_hist=True)
qx = np.linspace(np.min(q),np.max(q),1000)
qy = np.exp(tfd.Chi2(df=5).log_prob(qx))
sns.lineplot(qx,qy)
plt.show()
输出不是 DOF=5 的卡方分布。确实检验统计量经常有负值，这意味着优化后的结果往往比零假设拟合得更差，这应该是不可能的。
编辑 2：
这是对 "monster" 解决方案的尝试，其中我同时最小化了每个伪数据实现的不同输入变量的巨型网络。这感觉更像是 TensorFlow 可能擅长做的事情，尽管我觉得一旦我处理大量伪数据，我就会 运行 内存不足。尽管如此，我还是可以遍历成批的伪数据。
poises0 = [tfp.distributions.Poisson(rate = 10) for i in range(5)]
joint0 = tfd.JointDistributionSequential(poises0)
N = int(5e3)
samples0 = joint0.sample(N)
mus = [tf.Variable(10*np.ones(N, dtype='float32'), name='mu{0}'.format(i)) for i in range(5)]
poises_free = [tfp.distributions.Poisson(rate = mus[i]) for i in range(5)]
joint_free = tfd.JointDistributionSequential(poises_free)
qM = -2*(joint_free.log_prob(samples0))
@tf.function
def loss():
return tf.math.reduce_sum(qM,axis=0)
# Minimise
opt = tf.optimizers.SGD(0.1).minimize(loss,var_list=mus)
print("parameters:", mus)
print("loss:", loss())
q0 =-2*joint0.log_prob(samples0)
print("q0:", q0)
print("qM:", qM)
q = q0 - qM
fig = plt.figure()
ax = fig.add_subplot(111)
sns.distplot(q, kde=False, ax=ax, norm_hist=True)
qx = np.linspace(np.min(q),np.max(q),1000)
qy = np.exp(tfd.Chi2(df=5).log_prob(qx))
sns.lineplot(qx,qy)
plt.show()
不幸的是我现在得到错误：
Traceback (most recent call last):
File "testing3.py", line 35, in <module>
opt = tf.optimizers.SGD(0.1).minimize(loss,var_list=mus)
File "/home/farmer/anaconda3/envs/general/lib/python3.6/site-packages/tensorflow/python/keras/optimizer_v2/optimizer_v2.py", line 298, in minimize
return self.apply_gradients(grads_and_vars, name=name)
File "/home/farmer/anaconda3/envs/general/lib/python3.6/site-packages/tensorflow/python/keras/optimizer_v2/optimizer_v2.py", line 396, in apply_gradients
grads_and_vars = _filter_grads(grads_and_vars)
File "/home/farmer/anaconda3/envs/general/lib/python3.6/site-packages/tensorflow/python/keras/optimizer_v2/optimizer_v2.py", line 924, in _filter_grads
([v.name for _, v in grads_and_vars],))
ValueError: No gradients provided for any variable: ['mu0:0', 'mu1:0', 'mu2:0', 'mu3:0', 'mu4:0'].
我想这是一种基本的错误。我想我只是不明白 TensorFlow 如何跟踪它需要计算的导数。如果我在损失函数内部而不是外部定义变量，似乎一切正常，但我需要在外部使用它们以便稍后访问它们的值。所以我想我在这里有点不明白。
好的，这就是我的想法。我缺少的关键是：
- 将输入变量定义为巨大的张量，这样所有的最小化都可以同时发生。
- 一次为所有最小化构建一个组合损失函数
- 在损失函数定义中构建用于损失计算的中间变量，以便 TensorFlow 可以跟踪梯度（我认为
minimize 函数将损失函数包装在梯度带或类似的东西中）。
- 将损失函数定义为 class 的一部分，以便可以存储中间变量。
minimize只做了一步最小化，所以我们需要循环很多次直到它根据某种标准收敛。- 由于泊松分布均值小于零的无效性，我 运行 陷入了一些 NaN。所以我需要为输入变量添加一个约束。
有了这个，我现在可以在我的笔记本电脑上在大约 10 秒内完成相当于一百万次的最小化，这太棒了！
import tensorflow as tf
import tensorflow_probability as tfp
from tensorflow_probability import distributions as tfd
import seaborn as sns
import numpy as np
import matplotlib.pyplot as plt
# Bunch of independent Poisson distributions that we want to combine
poises0 = [tfd.Poisson(rate = 10) for i in range(5)]
# Construct joint distributions
joint0 = tfd.JointDistributionSequential(poises0)
N = int(1e6)
samples0 = joint0.sample(N)
class Model(object):
def __init__(self):
self.mus = [tf.Variable(10*np.ones(N, dtype='float32'), name='mu{0}'.format(i),
constraint=lambda x: tf.clip_by_value(x, 0.000001, np.infty)) for i in range(5)]
def loss(self):
poises_free = [tfd.Poisson(rate = self.mus[i]) for i in range(5)]
joint_free = tfd.JointDistributionSequential(poises_free)
# Construct (half of) test statistic
self.qM = -2*(joint_free.log_prob(samples0))
self.last_loss = tf.math.reduce_sum(self.qM,axis=0)
return self.last_loss
model = Model()
# Minimise
tol = 0.01 * N
delta_loss = 1e99
prev_loss = 1e99
i = 0
print("tol:", tol)
while delta_loss > tol:
opt = tf.optimizers.SGD(0.1).minimize(model.loss,var_list=model.mus)
delta_loss = np.abs(prev_loss - model.last_loss)
print("i:", i," delta_loss:", delta_loss)
i+=1
prev_loss = model.last_loss
q0 =-2*joint0.log_prob(samples0)
q = q0 - model.qM
print("parameters:", model.mus)
print("loss:", model.last_loss)
print("q0:", q0)
print("qM:", model.qM)
fig = plt.figure()
ax = fig.add_subplot(111)
sns.distplot(q, kde=False, ax=ax, norm_hist=True)
qx = np.linspace(np.min(q),np.max(q),1000)
qy = np.exp(tfd.Chi2(df=5).log_prob(qx))
sns.lineplot(qx,qy)
plt.show()
输出：
tol: 10000.0
i: 0 delta_loss: inf
i: 1 delta_loss: 197840.0
i: 2 delta_loss: 189366.0
i: 3 delta_loss: 181456.0
i: 4 delta_loss: 174040.0
i: 5 delta_loss: 167042.0
i: 6 delta_loss: 160448.0
i: 7 delta_loss: 154216.0
i: 8 delta_loss: 148310.0
i: 9 delta_loss: 142696.0
i: 10 delta_loss: 137352.0
i: 11 delta_loss: 132268.0
i: 12 delta_loss: 127404.0
...
i: 69 delta_loss: 11894.0
i: 70 delta_loss: 11344.0
i: 71 delta_loss: 10824.0
i: 72 delta_loss: 10318.0
i: 73 delta_loss: 9860.0
parameters: [<tf.Variable 'mu0:0' shape=(1000000,) dtype=float32, numpy=
array([ 6.5849004, 14.81182 , 7.506216 , ..., 10. , 11.491933 ,
10.760278 ], dtype=float32)>, <tf.Variable 'mu1:0' shape=(1000000,) dtype=float32, numpy=
array([12.881036, 7.506216, 12.881036, ..., 7.506216, 14.186232,
10.760278], dtype=float32)>, <tf.Variable 'mu2:0' shape=(1000000,) dtype=float32, numpy=
array([16.01586 , 8.378036 , 12.198007 , ..., 6.5849004, 12.198007 ,
8.378036 ], dtype=float32)>, <tf.Variable 'mu3:0' shape=(1000000,) dtype=float32, numpy=
array([10. , 7.506216, 12.198007, ..., 9.207426, 10.760278,
11.491933], dtype=float32)>, <tf.Variable 'mu4:0' shape=(1000000,) dtype=float32, numpy=
array([ 8.378036 , 14.81182 , 10. , ..., 6.5849004, 12.198007 ,
10.760278 ], dtype=float32)>]
loss: tf.Tensor(20760090.0, shape=(), dtype=float32)
q0: tf.Tensor([31.144037 31.440613 25.355555 ... 24.183338 27.195362 22.123463], shape=(1000000,), dtype=float32)
qM: tf.Tensor([21.74377 21.64162 21.526024 ... 19.488544 22.40428 21.08519 ], shape=(1000000,), dtype=float32)
结果现在卡方 DOF=5！或者至少非常接近。
我的最终目标是模拟似然比检验统计，然而，我遇到的核心问题是我不明白如何让 TensorFlow 2 对不同的数据输入执行许多优化。这是我的尝试，希望它能让您了解我正在尝试的事情：
import tensorflow as tf
import tensorflow_probability as tfp
from tensorflow_probability import distributions as tfd
import numpy as np
# Bunch of independent Poisson distributions that we want to combine
poises0 = [tfp.distributions.Poisson(rate = 10) for i in range(5)]
# Construct joint distributions
joint0 = tfd.JointDistributionSequential(poises0)
# Generate samples
N = int(1e3)
samples0 = joint0.sample(N)
# Now we need the same distributions but with floating parameters,
# and need to define the function to be minimised
mus = [tf.Variable(np.random.randn(), name='mu{0}'.format(i)) for i in range(5)]
#@tf.function
def loss():
poises_free = [tfp.distributions.Poisson(rate = mus[i]) for i in range(5)]
joint_free = tfd.JointDistributionSequential(poises_free)
# Construct (half of) test statistic
return -2*(joint_free.log_prob(samples0))
# Minimise (for all samples? Apparently not?)
opt = tf.optimizers.SGD(0.1).minimize(loss,var_list=mus)
print(mus)
print(loss())
print(opt)
quit()
输出：
[<tf.Variable 'mu0:0' shape=() dtype=float32, numpy=53387.016>, <tf.Variable 'mu1:0' shape=() dtype=float32, numpy=2540.568>, <tf.Variable 'mu2:0' shape=() dtype=float32, numpy=-5136.6226>, <tf.Variable 'mu3:0' shape=() dtype=float32, numpy=-3714.5227>, <tf.Variable 'mu4:0' shape=() dtype=float32, numpy=1062.9396>]
tf.Tensor(
[nan nan nan nan ... nan nan nan], shape=(1000,), dtype=float32)
<tf.Variable 'UnreadVariable' shape=() dtype=int64, numpy=1>
最后我想计算检验统计量
q = -2*joint0.log_prob(samples0) - loss()
并表明它服从具有 5 个自由度的卡方分布。
我是 TensorFlow 的新手，所以也许我做的完全错了，但我希望你明白我想要什么。
编辑：
所以我又玩了一会儿，我想 TensorFlow 根本不会像我假设的那样对输入张量进行并行优化。或者也许可以，但我需要进行不同的设置，即可能一次给它一个输入参数的张量和一个巨大的联合损失函数来进行所有最小化？
我也试过用一个简单的循环来做事情，看看会发生什么。正如预测的那样，它慢得可怜，但我什至没有得到正确的答案：
poises0 = [tfp.distributions.Poisson(rate = 10) for i in range(5)]
joint0 = tfd.JointDistributionSequential(poises0)
N = int(5e2)
samples0 = joint0.sample(N)
mus = [tf.Variable(10., name='mu{0}'.format(i)) for i in range(5)]
#@tf.function
def loss(xi):
def loss_inner():
poises_free = [tfp.distributions.Poisson(rate = mus[i]) for i in range(5)]
joint_free = tfd.JointDistributionSequential(poises_free)
# Construct (half of) test statistic
return -2*(joint_free.log_prob(xi))
return loss_inner
# Minimise
# I think I have to loop over the samples... bit lame. Can perhaps parallelise though.
q = []
for i in range(N):
xi = [x[i] for x in samples0]
opt = tf.optimizers.SGD(0.1).minimize(loss=loss(xi),var_list=mus)
q += [-2*joint0.log_prob(xi) - loss(xi)()]
fig = plt.figure()
ax = fig.add_subplot(111)
sns.distplot(q, kde=False, ax=ax, norm_hist=True)
qx = np.linspace(np.min(q),np.max(q),1000)
qy = np.exp(tfd.Chi2(df=5).log_prob(qx))
sns.lineplot(qx,qy)
plt.show()
输出不是 DOF=5 的卡方分布。确实检验统计量经常有负值，这意味着优化后的结果往往比零假设拟合得更差，这应该是不可能的。
编辑 2：
这是对 "monster" 解决方案的尝试，其中我同时最小化了每个伪数据实现的不同输入变量的巨型网络。这感觉更像是 TensorFlow 可能擅长做的事情，尽管我觉得一旦我处理大量伪数据，我就会 运行 内存不足。尽管如此，我还是可以遍历成批的伪数据。
poises0 = [tfp.distributions.Poisson(rate = 10) for i in range(5)]
joint0 = tfd.JointDistributionSequential(poises0)
N = int(5e3)
samples0 = joint0.sample(N)
mus = [tf.Variable(10*np.ones(N, dtype='float32'), name='mu{0}'.format(i)) for i in range(5)]
poises_free = [tfp.distributions.Poisson(rate = mus[i]) for i in range(5)]
joint_free = tfd.JointDistributionSequential(poises_free)
qM = -2*(joint_free.log_prob(samples0))
@tf.function
def loss():
return tf.math.reduce_sum(qM,axis=0)
# Minimise
opt = tf.optimizers.SGD(0.1).minimize(loss,var_list=mus)
print("parameters:", mus)
print("loss:", loss())
q0 =-2*joint0.log_prob(samples0)
print("q0:", q0)
print("qM:", qM)
q = q0 - qM
fig = plt.figure()
ax = fig.add_subplot(111)
sns.distplot(q, kde=False, ax=ax, norm_hist=True)
qx = np.linspace(np.min(q),np.max(q),1000)
qy = np.exp(tfd.Chi2(df=5).log_prob(qx))
sns.lineplot(qx,qy)
plt.show()
不幸的是我现在得到错误：
Traceback (most recent call last):
File "testing3.py", line 35, in <module>
opt = tf.optimizers.SGD(0.1).minimize(loss,var_list=mus)
File "/home/farmer/anaconda3/envs/general/lib/python3.6/site-packages/tensorflow/python/keras/optimizer_v2/optimizer_v2.py", line 298, in minimize
return self.apply_gradients(grads_and_vars, name=name)
File "/home/farmer/anaconda3/envs/general/lib/python3.6/site-packages/tensorflow/python/keras/optimizer_v2/optimizer_v2.py", line 396, in apply_gradients
grads_and_vars = _filter_grads(grads_and_vars)
File "/home/farmer/anaconda3/envs/general/lib/python3.6/site-packages/tensorflow/python/keras/optimizer_v2/optimizer_v2.py", line 924, in _filter_grads
([v.name for _, v in grads_and_vars],))
ValueError: No gradients provided for any variable: ['mu0:0', 'mu1:0', 'mu2:0', 'mu3:0', 'mu4:0'].
我想这是一种基本的错误。我想我只是不明白 TensorFlow 如何跟踪它需要计算的导数。如果我在损失函数内部而不是外部定义变量，似乎一切正常，但我需要在外部使用它们以便稍后访问它们的值。所以我想我在这里有点不明白。
好的，这就是我的想法。我缺少的关键是：
- 将输入变量定义为巨大的张量，这样所有的最小化都可以同时发生。
- 一次为所有最小化构建一个组合损失函数
- 在损失函数定义中构建用于损失计算的中间变量，以便 TensorFlow 可以跟踪梯度（我认为
minimize函数将损失函数包装在梯度带或类似的东西中）。 - 将损失函数定义为 class 的一部分，以便可以存储中间变量。
minimize只做了一步最小化，所以我们需要循环很多次直到它根据某种标准收敛。- 由于泊松分布均值小于零的无效性，我 运行 陷入了一些 NaN。所以我需要为输入变量添加一个约束。
有了这个，我现在可以在我的笔记本电脑上在大约 10 秒内完成相当于一百万次的最小化，这太棒了！
import tensorflow as tf
import tensorflow_probability as tfp
from tensorflow_probability import distributions as tfd
import seaborn as sns
import numpy as np
import matplotlib.pyplot as plt
# Bunch of independent Poisson distributions that we want to combine
poises0 = [tfd.Poisson(rate = 10) for i in range(5)]
# Construct joint distributions
joint0 = tfd.JointDistributionSequential(poises0)
N = int(1e6)
samples0 = joint0.sample(N)
class Model(object):
def __init__(self):
self.mus = [tf.Variable(10*np.ones(N, dtype='float32'), name='mu{0}'.format(i),
constraint=lambda x: tf.clip_by_value(x, 0.000001, np.infty)) for i in range(5)]
def loss(self):
poises_free = [tfd.Poisson(rate = self.mus[i]) for i in range(5)]
joint_free = tfd.JointDistributionSequential(poises_free)
# Construct (half of) test statistic
self.qM = -2*(joint_free.log_prob(samples0))
self.last_loss = tf.math.reduce_sum(self.qM,axis=0)
return self.last_loss
model = Model()
# Minimise
tol = 0.01 * N
delta_loss = 1e99
prev_loss = 1e99
i = 0
print("tol:", tol)
while delta_loss > tol:
opt = tf.optimizers.SGD(0.1).minimize(model.loss,var_list=model.mus)
delta_loss = np.abs(prev_loss - model.last_loss)
print("i:", i," delta_loss:", delta_loss)
i+=1
prev_loss = model.last_loss
q0 =-2*joint0.log_prob(samples0)
q = q0 - model.qM
print("parameters:", model.mus)
print("loss:", model.last_loss)
print("q0:", q0)
print("qM:", model.qM)
fig = plt.figure()
ax = fig.add_subplot(111)
sns.distplot(q, kde=False, ax=ax, norm_hist=True)
qx = np.linspace(np.min(q),np.max(q),1000)
qy = np.exp(tfd.Chi2(df=5).log_prob(qx))
sns.lineplot(qx,qy)
plt.show()
输出：
tol: 10000.0
i: 0 delta_loss: inf
i: 1 delta_loss: 197840.0
i: 2 delta_loss: 189366.0
i: 3 delta_loss: 181456.0
i: 4 delta_loss: 174040.0
i: 5 delta_loss: 167042.0
i: 6 delta_loss: 160448.0
i: 7 delta_loss: 154216.0
i: 8 delta_loss: 148310.0
i: 9 delta_loss: 142696.0
i: 10 delta_loss: 137352.0
i: 11 delta_loss: 132268.0
i: 12 delta_loss: 127404.0
...
i: 69 delta_loss: 11894.0
i: 70 delta_loss: 11344.0
i: 71 delta_loss: 10824.0
i: 72 delta_loss: 10318.0
i: 73 delta_loss: 9860.0
parameters: [<tf.Variable 'mu0:0' shape=(1000000,) dtype=float32, numpy=
array([ 6.5849004, 14.81182 , 7.506216 , ..., 10. , 11.491933 ,
10.760278 ], dtype=float32)>, <tf.Variable 'mu1:0' shape=(1000000,) dtype=float32, numpy=
array([12.881036, 7.506216, 12.881036, ..., 7.506216, 14.186232,
10.760278], dtype=float32)>, <tf.Variable 'mu2:0' shape=(1000000,) dtype=float32, numpy=
array([16.01586 , 8.378036 , 12.198007 , ..., 6.5849004, 12.198007 ,
8.378036 ], dtype=float32)>, <tf.Variable 'mu3:0' shape=(1000000,) dtype=float32, numpy=
array([10. , 7.506216, 12.198007, ..., 9.207426, 10.760278,
11.491933], dtype=float32)>, <tf.Variable 'mu4:0' shape=(1000000,) dtype=float32, numpy=
array([ 8.378036 , 14.81182 , 10. , ..., 6.5849004, 12.198007 ,
10.760278 ], dtype=float32)>]
loss: tf.Tensor(20760090.0, shape=(), dtype=float32)
q0: tf.Tensor([31.144037 31.440613 25.355555 ... 24.183338 27.195362 22.123463], shape=(1000000,), dtype=float32)
qM: tf.Tensor([21.74377 21.64162 21.526024 ... 19.488544 22.40428 21.08519 ], shape=(1000000,), dtype=float32)
结果现在卡方 DOF=5！或者至少非常接近。