v__0.157

mec/mx.py

@@ -24,7 +24,7 @@ def beta_gmm(Y_i,X_i_k,Z_i_l,W_l_l ):
     return beta_k,objval
 
 
-def pi_inv(pi_t_y,theLambda_k_k,epsilon_t_i_k, xi_k_y,maxit = 100000, reltol=1E-8):
+def pi_inv(pi_t_y,theLambda_k_k,epsilon_t_i_k, xi_k_y,maxit = 100000, reltol=1E-8, require_grad =False):
     (K ,Y ) = xi_k_y.shape
     (T,I,_) = epsilon_t_i_k.shape
     n_t_i = np.ones((T,1)) @ np.ones((1,I)) / I
@@ -40,14 +40,17 @@ def pi_inv(pi_t_y,theLambda_k_k,epsilon_t_i_k, xi_k_y,maxit = 100000, reltol=1E-
         else:
             U_t_y = Up_t_y
 
-    pi_t_i_y = np.concatenate( [ np.exp(U_t_y[:,None,:]+ varepsilon_t_i_y - u_t_i[:,:,None] ), 
-                                np.exp( - u_t_i)[:,:,None]],axis=2)
-
-    Sigma = sp.kron(sp.eye(T),sp.bmat([[sp.kron( sp.eye(I),      np.ones((1,Y+1)))            ],
-                                   [sp.kron(np.ones((1,I)),  sp.diags([1],shape=(Y,Y+1)))]]) )
-    Deltapi = sp.diags(pi_t_i_y.flatten())
-    proj = sp.kron(sp.eye(T),sp.kron( sp.eye(I), sp.diags([1],shape=(Y+1,Y)).toarray()) )
-    A = (Sigma @ Deltapi @ Sigma.T).tocsc()
-    B = (Sigma @ Deltapi @ proj @ sp.kron( epsilon_t_i_k.reshape((-1,K)) , xi_k_y.T )).tocsc()
-    dUdLambda_t_y_k_k = - sp.linalg.spsolve(A,B).toarray().reshape((T,I+Y,K,K))[:,-Y:,:,:]
+    if require_grad:
+        pi_t_i_y = np.concatenate( [ np.exp(U_t_y[:,None,:]+ varepsilon_t_i_y - u_t_i[:,:,None] ), 
+                                    np.exp( - u_t_i)[:,:,None]],axis=2)
+
+        Sigma = sp.kron(sp.eye(T),sp.bmat([[sp.kron( sp.eye(I),      np.ones((1,Y+1)))            ],
+                                    [sp.kron(np.ones((1,I)),  sp.diags([1],shape=(Y,Y+1)))]]) )
+        Deltapi = sp.diags(pi_t_i_y.flatten())
+        proj = sp.kron(sp.eye(T),sp.kron( sp.eye(I), sp.diags([1],shape=(Y+1,Y)).toarray()) )
+        A = (Sigma @ Deltapi @ Sigma.T).tocsc()
+        B = (Sigma @ Deltapi @ proj @ sp.kron( epsilon_t_i_k.reshape((-1,K)) , xi_k_y.T )).tocsc()
+        dUdLambda_t_y_k_k = - sp.linalg.spsolve(A,B).toarray().reshape((T,I+Y,K,K))[:,-Y:,:,:]
+    else:
+        dUdLambda_t_y_k_k = None
     return(U_t_y, dUdLambda_t_y_k_k)

setup.py

@@ -2,7 +2,7 @@
 
 setup(
 	name="mec",
-	version="0.156",
+	version="0.157",
 	authors=["Alfred Galichon"],
 	author_email="[email protected]",
 	licence="",