+- HP Forums (https://www.hpmuseum.org/forum)
+-- Forum: HP Calculators (and very old HP Computers) (/forum-3.html)
+--- Forum: HP Prime (/forum-5.html)
+--- Thread: Bo-Cox for multivariate normalization (/thread-22240.html)
Bo-Cox for multivariate normalization - robmio - 08-29-2024 09:36 AM
Good morning, everyone. I must solve the following “optimization” problem: I have a matrix of n=6 rows and p=2 columns:
M9:=[[1,4],[2,8],[7,12],[8,16],[15,21],[19,3]]
Assuming that the following matrix is not distributed according to the multivariate normal, there is the Box-Cox method to normalize the matrix. Box and Cox (1964) present the family of transformations defined by:
Z = M9(i,j)^L(j) if L<>0
Z = ln(M9(i,j)) if L(j) = 0
where, in this case, L(j):=[L1,L2]. The precise criterion for transforming to multivariate normality is to find the collection of transformations L(j):=[L1,L2] such that the function:
F(L) = (-n/2)*ln(ABS(S(L)))+Σ((L(j)-1)*Σ(ln(M9(i,j)),i,1,n),j,1,p)
is maximized. It is therefore a question of maximizing F(L).
Can anyone help me? Thanks so much
RE: Bo-Cox for multivariate normalization - robmio - 08-29-2024 10:39 AM
The proposed problem, for a single column matrix, has been solved with the Newton Raphson method. In the program below, "gij" is the gradient, "Hij" the Hessian matrix.
I find it difficult to extend the program to matrices with more than one column (multivariate matrices). Is there perhaps a command in HP Prime to optimize the function described in the first post?
Code:
#cas
BoxCox_Newton_Cas(r,λo,iter):=
BEGIN
LOCAL y, z, u, v, x, Ii, Hh, sσ2;
LOCAL dmsr, j, gij, Hij;
y:=mat2list(r);
dmsr:=rowDim(r);
x:=MAKEMAT(1,dmsr,1);
Hh:=x*inv(TRN(x)*x)*TRN(x);
Ii:=IDENMAT(dmsr);
FOR j FROM 1 TO iter DO
IF λo≠0 THEN
z:=(y.^λo-1)./λo;
ELSE
z:=LN(y);
END;
IF λo≠0 THEN
u:=(λo*y.^λo*LN(y)-y.^λo+1)./λo.^2;
ELSE
u:=LN(y).^2/2;
END;
IF λo≠0 THEN
v:=(2-((λo*LN(y)).^2-2*λo*LN(y)+2)*
y.^λo)./λo.^3;
ELSE
v:=−LN(y).^3/3;
END;
z:=list2mat(z,1);
sσ2:=(TRN(z)*(Ii-Hh)*z)/dmsr;
u:=list2mat(u,1);
gij:=−(TRN(u)*(Ii-Hh)*z)*inv(sσ2)+
ΣLIST(LN(y));
v:=list2mat(v,1);
Hij:=(TRN(v)*(Ii-Hh)*z-TRN(u)*(Ii-Hh)*
u)*inv(sσ2)+2/dmsr*((TRN(u)*(Ii-Hh)*
z)*inv(sσ2)).^2;
λo:=λo-inv(Hij)*gij;
λo:=GET(GET(λo,1),1);
END;
RETURN {λo,gij};
END;
#end