FYI, I suspect this doesn't do what you think it does. np.rank is the same as array.ndims IIRC. You want np.linalg.matrix_rank. Been a pain point for a while. We used to have our own version of rank that we used. 

This fixed the problem you were seeing with f_pvalue for the between model, but now I need to check the other degrees of freedom to see what else is going on.