def sum_outer_product_balanced(x, n_groups):
    xrs = x.reshape(-1, n_groups, order='F')
    return np.dot(xrs, xrs.T)  #should be (nobs_i, nobs_i)

    #x.reshape(n_groups, nobs_i,  k_vars) #, order='F')
    #... ? this is getting 3-dimensional  dot, tensordot?
    #needs (n_groups, k_vars, k_vars) array with sum over groups
    #NOT
    #I only need this for x is 1d, i.e. residual
sum outerproduct dot(x_i, x_i.T) over individuals
dataset/ETHPy150Open statsmodels/statsmodels/statsmodels/sandbox/panel/panel_short.py/sum_outer_product_balanced