function q = GMM(bgmm) global n dc dc1 dc2 dc3 dc4 dy; m1 = dc1.*(dc - bgmm*dy); m2 = dc2.*(dc - bgmm*dy); m3 = dc3.*(dc - bgmm*dy); m4 = dc4.*(dc - bgmm*dy); m = [m1,m2,m3,m4]; W = (1/n)*inv(m'*m); mbar = mean(m)'; q = mbar'*W*mbar; end