吴恩达的机器学习和深度学习都讲得很好,得亏B站大佬们资源丰富和无私奉献,视频、笔记和课后测验都有,我这也是给自己留个门。课后测验题参考两位大佬阿宽的博客–超详细课后测验和大大鹏I6的B站Python视频资料汇总。
前情提要
- 二分类问题:二分类算法最简单的例子就是猫的识别问题,给定一张图片作为输入,识别图片是否为猫(1和0作为输出标签表示)。一张图片在计算机中使用三个矩阵(红、绿、蓝三种颜色)作为图片的特征信息来保存,为了简化可以将三个矩阵连接成一个特征向量,这样每个样本可以用一个输入特征向量x和输出标签y表示,如此整个训练集可以表示成更简单的训练集矩阵(每一列都是一个样本的特征向量)。
- 逻辑回归算法:适用于二分类问题。就是针对给定一个样本的输入特征向量,如何去表示这个样本输出标签的估计值,如果直接用线性函数
y
^
=
w
T
x
+
b
,
\hat y = w^{T}x+b,
y^=wTx+b,表示,容易出现该估计值
y
^
\hat y
y^ 不是介于0和1之间(因为用0到1之间的数代替
y
^
\hat y
y^表示实际值y等于1的概率),所有在线性函数的基础上增加一个激活函数sigmoid函数,即令
z
=
w
T
x
+
b
,
z= w^{T}x+b,
z=wTx+b, 而最终的输出标签
y
^
=
σ
(
z
)
\hat y=\sigma(z)
y^=σ(z)代替,如此便可以使得满足结果在0和1之间。
实现逻辑回归时的主要工作就是训练机器学习参数 w 和 b 这样才能使得 y ^ \hat y y^成为对 y=1 这一情况的很好估计。另外,可以更简化该表达式,令参数 θ \theta θ 表示参数w和b的连接向量,即 y ^ = σ ( θ T x ) \hat y= \sigma(\theta^{T}x) y^=σ(θTx)。 - 逻辑回归的代价函数:前面我们得到逻辑回归是线性函数和激活函数的复合,那么对于线性函数中的参数
w
w
w 和
b
b
b 就需要训练得到最佳的值,那么代价函数就是整个训练集的训练结果与实际结果之间的误差,是算法的在全部训练数据上的表现,所以通过找到最小的总代价
J
J
J来确定最佳的参数
w
w
w和
b
b
b。
(1) 损失函数——又叫误差函数,用于衡量算法的运行情况,预测输出值和实际值的接近程度,一般用两者的平方差或其一半。但在逻辑回归模型中我们定义另外一个损失函数 L o s t F u n c t i o n Lost \ Function Lost Function: L ( y ^ , y ) = − y l o g ( y ^ ) − ( 1 − y ) l o g ( 1 − y ^ ) , L(\hat y,y)=-ylog(\hat y)-(1-y)log(1-\hat y), L(y^,y)=−ylog(y^)−(1−y)log(1−y^), 这样定义损失函数的最简单的原因就是当 y = 1 y=1 y=1时, y ^ \hat y y^越大则损失函数 L L L越小,同理当 y = 0 y=0 y=0时, y ^ \hat y y^越小则损失函数 L L L越小,即保证 y ^ \hat y y^和实际值 y y y是保持同向变化的;
(2) 代价函数 J J J——损失函数是单个训练样本中定义的误差,而代价函数就是对整个训练集m个样本的损失函数求均值: J ( w , b ) = ∑ i = 1 m L ( y ^ ( i ) , y ( i ) ) , J(w,b)=\sum_{i=1}^{m} L(\hat y^{(i)},y^{(i)}), J(w,b)=i=1∑mL(y^(i),y(i)),其中损失函数 L L L的定义同上,而 y ( i ) y^{(i)} y(i)表示第 i i i个样本的实际结果。 - 梯度下降法:梯度下降法就是找到使得代价函数
J
J
J取最小值的对应的
w
w
w和
b
b
b,基本思想就是朝着梯度下降最快的方向迭代直至找到最佳结果。首先随机初始化参数
w
w
w和
b
b
b(因为
J
J
J是凸函数,所以初始化的位置在哪都可以到达大致同一个位置),然后朝梯度下降最快的方向(最陡的下坡方向不断迭代)直到走到(或接近)全局最优解的位置,如下图所示。
而梯度下降法中更需要注意的是参数 w w w和 b b b的迭代更新,即如逻辑回归中的代价函数 J ( w , b ) J(w,b) J(w,b)有两个参数,每次迭代过程中的参数的值变化: w : = w − α ∂ J ( w , b ) ∂ w , w:=w - \alpha \frac{\partial J(w,b)}{\partial w}, w:=w−α∂w∂J(w,b), b : = b − α ∂ J ( w , b ) ∂ b b:=b - \alpha \frac{\partial J(w,b)}{\partial b} b:=b−α∂b∂J(w,b)利用偏导数来表示梯度,而 α \alpha α是学习率,表示步长。 - 计算图:其实就是一个神经网络的计算流程图,包含正向从左到右计算输出结果的过程和反向从右到左计算各个导数过程。正向过程就类似于将逐步运算用树的形式画出来,而反向传播过程就是计算导数,也即是后一个阶段中的节点结果随前一个阶段中的节点结果变化的快慢,通过计算图可以很清晰地看到其中的变化关系。具体的计算图就类似于如下这样的流程图。
- 逻辑回归中的梯度下降法:就是通过计算偏导数来实现逻辑回归中的梯度下降法。因为梯度下降法中,梯度的求解是主要部分,需用来更新参数,而在逻辑回归的代价函数中参数有
w
w
w和
b
b
b,因此梯度选择偏导数来求解。假设输入特征为
x
1
x_1
x1和
x
2
x_2
x2,那么对应就需要参数
w
=
[
w
1
,
w
2
]
w=[w_1,w_2]
w=[w1,w2]和参数
b
b
b,因此有中间结果
z
:
=
w
1
x
1
+
w
2
x
2
+
b
=
w
T
x
+
b
z:=w_1x_1+w_2x_2+b=w^Tx+b
z:=w1x1+w2x2+b=wTx+b (输入特征变量的线性组合),预测结果标签
y
^
=
σ
(
z
)
=
1
1
+
e
−
z
\hat y=\sigma(z)=\frac{1}{1+e^{-z}}
y^=σ(z)=1+e−z1 (利用
σ
(
x
)
\sigma(x)
σ(x)作为激活函数作用在线性组合的中间结果上),最后可以计算出代价函数:
J
(
w
,
b
)
=
1
m
∑
i
=
1
m
L
(
y
^
(
i
)
,
y
(
i
)
)
=
1
m
∑
i
=
1
m
−
y
(
i
)
l
o
g
y
^
(
i
)
−
(
1
−
y
(
i
)
)
l
o
g
(
1
−
y
^
(
i
)
)
.
J(w,b)=\frac{1}{m}\sum^m_{i=1}L(\hat y^{(i)},y^{(i)})=\frac{1}{m}\sum^m_{i=1}-y^{(i)}log\hat y^{(i)}-(1-y^{(i)})log(1-\hat y^{(i)}).
J(w,b)=m1i=1∑mL(y^(i),y(i))=m1i=1∑m−y(i)logy^(i)−(1−y(i))log(1−y^(i)).如果是单个样本的情况,其代价函数
L
(
a
,
y
)
=
−
(
y
l
o
g
(
a
)
+
(
1
−
y
)
l
o
g
(
1
−
a
)
)
,
L(a,y)=-(ylog(a)+(1-y)log(1-a)),
L(a,y)=−(ylog(a)+(1−y)log(1−a)),
a
a
a是逻辑回归的输出,
y
y
y是样本的标签值,
w
w
w和
b
b
b的修正量就是
w
:
=
w
−
α
∂
J
(
w
,
b
)
∂
w
,
b
:
=
b
−
α
∂
J
(
w
,
b
)
∂
b
。
w:=w-\alpha \frac{\partial J(w,b)}{\partial w},b:=b-\alpha \frac{\partial J(w,b)}{\partial b}。
w:=w−α∂w∂J(w,b),b:=b−α∂b∂J(w,b)。
依旧考虑单个样本的情况,根据计算图反向求解导数(梯度),我们的目的就是找到使得逻辑回归中代价函数的最小值 L ( a , y ) L(a,y) L(a,y),整个过程需要做的就是利用上述修正量公式更新参数 w w w和 b b b的值即可。首先预测结果 a a a是输入特征变量的线性组合,也是参数 w w w和 b b b的复合函数,那么代价函数 J J J对于参数 w w w和 b b b的梯度(导数)求解需要先求解对于 a a a的导数,即 d J ( w , b ) d a = d L ( a , y ) d a = − y a + 1 − y 1 − a \frac{dJ(w,b)}{da}=\frac{dL(a,y)}{da}=-\frac{y}{a}+\frac{1-y}{1-a} dadJ(w,b)=dadL(a,y)=−ay+1−a1−y,而 a = σ ( z ) a=\sigma(z) a=σ(z)可得 d a d z = e z ( e z + 1 ) 2 = a ( 1 − a ) \frac{da}{dz}=\frac{e^z}{(e^z+1)^2}=a(1-a) dzda=(ez+1)2ez=a(1−a),故 d J ( w , b ) d z = d L ( a , y ) d z = d L ( a , y ) d a ⋅ d a d z = ( − y a + 1 − y 1 − a ) ⋅ ( a ( 1 − a ) ) = a − y , \frac{dJ(w,b)}{dz}=\frac{dL(a,y)}{dz}=\frac{dL(a,y)}{da}·\frac{da}{dz}=(-\frac{y}{a}+\frac{1-y}{1-a})·(a(1-a))=a-y, dzdJ(w,b)=dzdL(a,y)=dadL(a,y)⋅dzda=(−ay+1−a1−y)⋅(a(1−a))=a−y,那么最后一步反向推导计算就是关于参数 w w w和 b b b的导数,因为 z z z是 w w w和 b b b的线性组合,所以易得 d z d w 1 = x 1 , d z d w 2 = x 2 , d z d w 2 = 1 。 \frac{dz}{dw_1}=x_1, \frac{dz}{dw_2}=x_2,\frac{dz}{dw_2}=1。 dw1dz=x1,dw2dz=x2,dw2dz=1。因此,可以得到代价函数对于参数 w w w和 b b b的导数(分别记为 d w 1 , d w 2 , d b dw_1, dw_2, db dw1,dw2,db)如下: d w 1 = d L ( a , y ) d w 1 = x 1 ⋅ d L ( a , y ) d z = x 1 ⋅ ( a − y ) , dw_1=\frac{dL(a,y)}{dw_1}=x_1·\frac{dL(a,y)}{dz}=x_1·(a-y), dw1=dw1dL(a,y)=x1⋅dzdL(a,y)=x1⋅(a−y), d w 2 = d L ( a , y ) d w 2 = x 2 ⋅ d L ( a , y ) d z = x 2 ⋅ ( a − y ) , dw_2=\frac{dL(a,y)}{dw_2}=x_2·\frac{dL(a,y)}{dz}=x_2·(a-y), dw2=dw2dL(a,y)=x2⋅dzdL(a,y)=x2⋅(a−y), d b = d L ( a , y ) d b = d L ( a , y ) d z = ( a − y ) , db=\frac{dL(a,y)}{db}=\frac{dL(a,y)}{dz}=(a-y), db=dbdL(a,y)=dzdL(a,y)=(a−y),代入到参数的修正更新公式中即可。至此,单个样本两个输入特征变量的情况就是这样计算的,如果是若干m个样本(仍是两个输入特征变量情况),那么计算结果的公式为: d w 1 = 1 m ∑ i = 1 m x 1 ( i ) ⋅ d L ( a ( i ) , y ( i ) ) d z = 1 m ∑ i = 1 m x 1 ( i ) ⋅ ( a ( i ) − y ( i ) ) , dw_1=\frac{1}{m}\sum_{i=1}^mx_1^{(i)}·\frac{dL(a^{(i)},y^{(i)})}{dz}=\frac{1}{m}\sum_{i=1}^mx_1^{(i)}·(a^{(i)}-y^{(i)}), dw1=m1i=1∑mx1(i)⋅dzdL(a(i),y(i))=m1i=1∑mx1(i)⋅(a(i)−y(i)),
d w 2 = 1 m ∑ i = 1 m x 2 ( i ) ⋅ d L ( a ( i ) , y ( i ) ) d z = 1 m ∑ i = 1 m x 2 ( i ) ⋅ ( a ( i ) − y ( i ) ) , dw_2=\frac{1}{m}\sum_{i=1}^mx_2^{(i)}·\frac{dL(a^{(i)},y^{(i)})}{dz}=\frac{1}{m}\sum_{i=1}^mx_2^{(i)}·(a^{(i)}-y^{(i)}), dw2=m1i=1∑mx2(i)⋅dzdL(a(i),y(i))=m1i=1∑mx2(i)⋅(a(i)−y(i)),
d b = 1 m ∑ i = 1 m d L ( a ( i ) , y ( i ) ) d z = 1 m ∑ i = 1 m ( a ( i ) − y ( i ) ) , db=\frac{1}{m}\sum_{i=1}^m\frac{dL(a^{(i)},y^{(i)})}{dz}=\frac{1}{m}\sum_{i=1}^m(a^{(i)}-y^{(i)}), db=m1i=1∑mdzdL(a(i),y(i))=m1i=1∑m(a(i)−y(i)),
再结合 w 1 : = w 1 − α ⋅ d w 1 , w 2 : = w 2 − α ⋅ d w 2 , b : = b − α ⋅ d b w_1:=w_1-\alpha·dw_1,w_2:=w_2-\alpha·dw_2,b:=b-\alpha·db w1:=w1−α⋅dw1,w2:=w2−α⋅dw2,b:=b−α⋅db 可得出两个参数的更新值,代入下一次迭代计算。
-m个样本的梯度下降:上面讨论的是2个输入特征变量的 m 个样本逻辑回归中的梯度下降计算公式,其中只是给出了2个输入特征变量 x 1 x_1 x1和 x 2 x_2 x2,所以对应在代价函数中参数 w w w是二维的,而在计算导数过程中需要相应地计算特征变量的个数次,所以对于上述结果每次迭代只需要执行两次(即参数 w w w的维数,或是 x x x的输入特征个数)。如果是m个样本的话,我们有以下的计算结果及公式:
d w k = 1 m ∑ i = 1 m x k ( i ) ⋅ d L ( a ( i ) , y ( i ) ) d z = 1 m ∑ i = 1 m x k ( i ) ⋅ ( a ( i ) − y ( i ) ) , dw_k=\frac{1}{m}\sum_{i=1}^mx_k^{(i)}·\frac{dL(a^{(i)},y^{(i)})}{dz}=\frac{1}{m}\sum_{i=1}^mx_k^{(i)}·(a^{(i)}-y^{(i)}), dwk=m1i=1∑mxk(i)⋅dzdL(a(i),y(i))=m1i=1∑mxk(i)⋅(a(i)−y(i)),
d b = 1 m ∑ i = 1 m d L ( a ( i ) , y ( i ) ) d z = 1 m ∑ i = 1 m ( a ( i ) − y ( i ) ) , db=\frac{1}{m}\sum_{i=1}^m\frac{dL(a^{(i)},y^{(i)})}{dz}=\frac{1}{m}\sum_{i=1}^m(a^{(i)}-y^{(i)}), db=m1i=1∑mdzdL(a(i),y(i))=m1i=1∑m(a(i)−y(i)),
根据前一条逻辑回归的GD方法中两个输入特征变量的计算公式以及上述m个样本的计算公式可以利用编码的方式实现,首先视频中给出了代码实现两个输入特征变量对应的逻辑回归的梯度下降方法。
自然更一般的对于m个样本、n个输入特征向量的代码实现,最重要的也就是反向传播过程中的导数(梯度)计算和对代价函数中两个参数 w w w和 b b b的更新,就是在两层循环,外循环是针对每个输入特征向量 x k x_k xk的参数 w k w_k wk的特征个数 n n n,然后内循环是计算每个参数 w k w_k wk在预测结果与真实结果的差值累加与求平均的过程。
首先初始化 J = 0 ; d w k = 0 ; d b = 0 J=0;dw_k=0;db=0 J=0;dwk=0;db=0,然后就是主体代码流程
如上我们得到的就只是一次迭代(一步梯度下降),所以该计算过程需要重复很多次才可以应用多次梯度下降而达到想要的效果,而且该过程需要编写两个for循环才能进行一次完整的梯度下降。因为 for 循环在深度学习中会使得算法很低效,那么就需要借用python算法的优势,就是向量化技术,将for循环转换成向量进行计算,从而减少循环层数提高深度学习的效果。for k=1 to n for i = 1 to m z(i) = w_k·x_k(i)+b; a(i) = σ(z(i)); J += -[y(i)log(a(i))+(1-y(i))log(1-a(i))]; dz(i) = a(i)-y(i); dw_k += x_k(i)dz(i); db += dz(i); J/= m; dw_k/= m; db/= m; w_k=w_k-α*dw_k b=b-α*db
- 向量化:向量化是基础的去除代码中for循环的方法,将循环转换为向量中的一维向量基,而在深度学习的训练中代码运行速度很重要,所以将循环转换为向量是一个很关键的技巧。例如在逻辑回归中计算输入特征变量的线性组合
z
=
w
T
⋅
x
+
b
z=w^T·x+b
z=wT⋅x+b,其中参数
w
w
w和输入特征
x
x
x都是列向量,而且参数
w
w
w是与输入特征$x相关的,如果使用非向量化的方法计算该线性组合则需要用for循环如下
上述过程的非向量化实现(python)对每个特征变量和对应的参数都需要通过循环来计算,而在python中可以通过向量化直接计算 w T ⋅ x w^T·x wT⋅x,利用x z=0 for i in range(n_x) z+=w[i]*x[i] z+=b
可能说这么多都是理论上的,我们利用Jupyter Notebook编写代码来比较向量化计算多维数组和for循环计算多维数组的时间差别。具体结果如下:z=np.dot(w,x)+b #其中w和x是列向量,w=[w_1,...,w_n], x=[x_1,...,x_n]
通过上述结果很明显能比较出向量化的计算时间比for循环的计算时间短300倍,这意味着向量化带来的计算效率的提升是很有影响的。import numpy as np #导入numpy库 a = np.array([1,2,3,4]) #创建一个数据a,感觉都不一样 print(a) [1]output:[1 2 3 4] import time #导入时间库 a = np.random.rand(1000000) b = np.random.rand(1000000) #通过round随机得到两个一百万维度的数组 tic = time.time() #现在测量一下当前时间 #向量化的版本 c = np.dot(a,b) toc = time.time() print(c) print("Vectorized version:" + str(1000*(toc-tic)) +"ms") #打印一下向量化的版本的时间 #继续增加非向量化的版本 c = 0 tic = time.time() for i in range(1000000): c += a[i]*b[i] toc = time.time() #循环结束记录下结束时间 print(c) print("For loop:" + str(1000*(toc-tic)) + "ms") #打印for循环的版本的时间 [2]output: 250222.34816827404 Vectorized version: 1.9948482513427734ms 250222.34816827823 For loop: 582.4785232543945ms
在jupyter notebook上面实现python计算只有CPU,而CPU和GPU都有并行化的指令,叫做SIMD指令----代表了一个单独指令多维数据,基础意义是,如果你使用了built-in函数,像np.function或者并不要求你实现循环的函数,它可以让python的充分利用并行化计算,在GPU和CPU上面计算,GPU更加擅长SIMD计算,但是CPU只是相比较而言较差一点。 - 向量化的更多例子:在上述的for循环和numpy库的向量化函数对比中,可以发现向量化的计算效率比for循环高不少,而在numpy库中有很多向量函数(都是内置函数),比如计算指数函数(exp)的
u=np.exp(v)
、计算对数函数(log)的u=np.log()
、计算数据绝对值np.abs()
、计算元素 v v v中的最大值np.maximum() or np.maximum(v, 0)
、v**2
代表获得元素 v v v每个值的平方、1/v
代表获取元素 v v v的倒数等等。将向量化利用到逻辑回归的梯度下降中,求导部分有两层循环,当输入特征变量有 n ( n ≥ 3 ) n(n≥3) n(n≥3)个时,需要循环计算 d w 1 、 d w 2 、 d w 3 dw_1、dw_2、dw_3 dw1、dw2、dw3等,那么可以利用向量化将里面对于特征的循环去掉。
- 向量化逻辑回归的梯度输出:将训练输入定义一个矩阵
X
X
X,将所有样本按列堆积起来形成一个
n
x
n_x
nx行
m
m
m列的矩阵(因为有m个样本,每个样本有
n
x
n_x
nx个元素),那么逻辑回归的过程可以表示如下
按上图所示,计算就是要计算每个 z ( i ) z^{(i)} z(i),其中可以利用向量化将所有的按列存放一个行向量,最终可以表示为 [ z ( 1 ) z ( 2 ) . . . z ( m ) ] = w T ⋅ X + [ b b . . . b ] = [ w T x ( 1 ) + b , w T x ( 2 ) + b , . . . , w T x ( m ) + b ] [z^{(1)} z^{(2)}... z^{(m)}]=w^{T}·X+[bb...b]=[w^{T}x^{(1)}+b,w^{T}x^{(2)}+b,...,w^{T}x^{(m)}+b] [z(1)z(2)...z(m)]=wT⋅X+[bb...b]=[wTx(1)+b,wTx(2)+b,...,wTx(m)+b]利用Python代码实现中间的向量计算使用numpy Z = n p . d o t ( w . T , X ) + b Z=np.dot(w.T, X)+b Z=np.dot(w.T,X)+b命令即可。这里b是一个实数,但是在向量计算时在python中会自动将这个实数扩展为一个 1 × m 1 \times m 1×m向量,这在python中称之为广播。在计算得到输入特征的组合后,还需要通过激活函数 σ ( z ) \sigma(z) σ(z),也是可以通过一行代码 A = [ a ( 1 ) a ( 2 ) . . . a ( m ) ] = σ ( Z ) A=[a^{(1)} a^{(2)}...a^{(m)}]=\sigma(Z) A=[a(1)a(2)...a(m)]=σ(Z)一次性计算得出所有 a a a。这就是在同一时间内完成所有m个训练样本的前向传播向量化计算。
上述我们已经有了向量化同时计算m个训练样本的前向传播计算,但是后向传播同时计算m个数据的梯度仍有一个对训练集的循环,也可以通过向量化计算利用Python中numpy库函数简化计算,对比如下
没有向量化
向量化
Z=w^T·X+b=np.dot(w.T,X)+b A=σ(Z) dZ=A-Y dw=(1/m)*X*d(z^T) db=(1/m)*np.sum(dZ) w:=w-α*dw b:=b-α*db
利用上述代码中的前五个公式完成了前向和后向传播,也实现了对所有训练样本进行预测和求导,再利用后两个公式对梯度下降更新参数,这样我们便利用向量化简化了整个逻辑回归的梯度下降算法。
课后测验–DL简介(第二周)
- 神经元节点的职责或功能
首先输入层的特征向量(矩阵),隐藏层的神经元节点先是计算特征向量的线性组合,再是利用激活函数计算激活,最后是输出标签。 - 假设 img 是一个(32,32,3)数组,具有3个颜色通道:红色、绿色和蓝色的32x32像素的图像。 如何将其重新转换为列向量?
利用reshape函数方法,将三个32×32的矩阵转换为32×32×3的列向量。
x=img.reshape((32 * 32 * 3, 1))
- 看一下下面的这两个随机数组“a”和“b”,计算数组c的维度
b(列向量)复制3次,以便它可以和A的每一列相加,所以:a = np.random.randn(2, 3) # a.shape = (2, 3) b = np.random.randn(2, 1) # b.shape = (2, 1) c = a + b
c.shape = (2, 3)
- 看一下下面的这两个随机数组“a”和“b”,请问数组“c”的维度是多少
运算符 “*” 说明了按元素乘法来相乘,但是元素乘法需要两个矩阵之间的维数相同,所以将报错无法计算。 - 假设你的每一个实例有
n
x
n_x
nx个输入特征,想一下在
X
=
[
x
(
1
)
,
x
(
2
)
…
x
(
m
)
]
X=[x^{(1)}, x^{(2)}…x^{(m)}]
X=[x(1),x(2)…x(m)]中,
X
X
X的维度是多少?
( n x , m ) (n_x, m) (nx,m) - 回想一下,
np.dot(a,b
在a和b上执行矩阵乘法,而`a * b’执行元素方式的乘法。 看一下下面的这两个随机数组“a”和“b”,请问c的维度是多少:a = np.random.randn(12288, 150) # a.shape = (12288, 150) b = np.random.randn(150, 45) # b.shape = (150, 45) c = np.dot(a, b)
c.shape = (12288, 45)
, 这是一个简单的矩阵乘法例子。 - 看一下下面的这个代码片段:
请问要怎么把它们向量化?# a.shape = (3,4) # b.shape = (4,1) for i in range(3): for j in range(4): c[i][j] = a[i][j] + b[j]
c = a + b.T