How to calculate partial derivatives of a function with respect to matrix elements?

Question

Python

INPUT_DIM = 2
OUT_DIM = 1
H_DIM = 2

W1 = np.random.randn(INPUT_DIM, H_DIM)
b1 = np.random.randn(H_DIM)
W2 = np.random.randn(H_DIM, OUT_DIM)
b2 = np.random.randn(OUT_DIM)

def train(data, all_trues):
  learn_rate = 0.1
  epochs = 1000
  for epoch in range(epochs):
    for x,y_true in zip(data,all_trues):
      t = x @ W1 + b1
      h = sigmoid(t)
      t2 = h @ W2 + b2
      o = sigmoid(t2)
      y_pred = o
      d_L_d_ypred = -2*(y_true - y_pred)
      d_h_d_w1 = diff(t,W1)*sigmoid(t) #можно ли как-нибудь посчитать все частные производные по элементам матрицы?