code: def chisquare(n_ij, weighted):
      
      if weighted:
          m_ij = n_ij / n_ij
  
          nan_mask = np.isnan(m_ij)
          m_ij[nan_mask] = 0.000001  # otherwise it breaks the chi-squared test
  
          w_ij = m_ij
          n_ij_col_sum = n_ij.sum(axis=1)
          n_ij_row_sum = n_ij.sum(axis=0)
          alpha, beta, eps = (1, 1, 1)
          while eps > 10e-6:
              alpha = alpha * np.vstack(n_ij_col_sum / m_ij.sum(axis=1))
              beta = n_ij_row_sum / (alpha * w_ij).sum(axis=0)
              eps = np.max(np.absolute(w_ij * alpha * beta - m_ij))
              m_ij = w_ij * alpha * beta
  
      else:
          m_ij = (np.vstack(n_ij.sum(axis=1)) * n_ij.sum(axis=0)) / n_ij.sum().astype(float)
  
      dof = (n_ij.shape[0] - 1) * (n_ij.shape[1] - 1)
      chi, p_val = stats.chisquare(n_ij, f_exp=m_ij, ddof=n_ij.size - 1 - dof, axis=None)
  
      return (chi, p_val, dof)
docs: Calculates the chisquare for a matrix of ind_v x dep_v
      for the unweighted and SPSS weighted case