HP Forums

Full Version: HP Prime for advance statistics
You're currently viewing a stripped down version of our content. View the full version with proper formatting.
Does anyone know of a one way anova and or multiple regression app that can be loaded to the HP Prime?
(09-04-2014 07:02 AM)pr011235 Wrote: [ -> ]Does anyone know of a one way anova and or multiple regression app that can be loaded to the HP Prime?

I am attaching a pdf file that contains several articles (from the HP Solve newsletter) for the HP 39gII that perform ADVANCED types of regression. These programs should work with the HP Prime. These programs show the power of the HP 39gII and the HP Prime in performing sophisticated regression calculations.

Namir

PS: Here is the code for about the dozen programs in the attached article, so you can copy and pasted them easily in the HP Prime emulator or the Connectivity Kit.

MLR2

Code:
EXPORT MLR2(MatX,VectY)
BEGIN
  LOCAL i,j,y,Rsqr;
  LOCAL lstDimX,NumRowsX;
  LOCAL YMean,Sum1,Sum2,Yhat,RegCoeff;
  // calculate the regression coefficients
  RegCoeff:=LSQ(MatX,VectY);
  // calculate the number of data points
  lstDimX:=SIZE(MatX);
  NumRowsX:=lstDimX(1);
  // calculate ymean
  Sum1:=0;
  FOR i FROM 1 TO NumRowsX DO
    y:=VectY(i,1);
    Sum1:=Sum1+y;
  END;
  YMean:=Sum1/NumRowsX;
  // calculate the coefficient of determination
  Sum1:=0;
  Sum2:=0;
  Yhat:=MatX*RegCoeff;
  FOR i FROM 1 TO NumRowsX DO
    Sum1:=Sum1+(Yhat(i,1)-YMean)^2;
    Sum2:=Sum2+(VectY(i,1)-YMean)^2;
  END;
  Rsqr:=Sum1/Sum2;
  // return the results
  RETURN {Rsqr,RegCoeff};
END;

PolyReg
Code:
EXPORT PolyReg(DSMat,SelXCol,SelYCol,Order)
BEGIN
  LOCAL i,j,x,y;
  LOCAL lstDimX,NumRows;
  LOCAL MatX, VectY, RegCoeff;
  LOCAL YMean,Sum1,Sum2,Yhat,Rsqr;
  lstDimX:=SIZE(DSMat);
  NumRows:=lstDimX(1);
  // initialize matrix and vector
  MatX:=MAKEMAT(1,NumRows,Order+1);
  VectY:=MAKEMAT(1,NumRows,1);
  Sum1:=0;
  // populate matrix X and vector y with data
  FOR i FROM 1 TO NumRows DO
    y:=DSMat(i,SelYCol);
    x:=DSMat(i,SelXCol);
    VectY(i,1):=y;
    Sum1:=Sum1+y;
    FOR j FROM 1 TO Order DO
      MatX(i,j+1):=x^j;
    END;
  END;
  // calculate ymean
  YMean:=Sum1/NumRows;
  // calculate regression coefficients
  RegCoeff:=LSQ(MatX,VectY);
  // calculate the sums of squares of
  // (yhat - ymean) and (y - ymean)
  Sum1:=0;
  Sum2:=0;
  FOR i FROM 1 TO NumRows DO
    y:=DSMat(i,SelYCol);
    x:=DSMat(i,SelXCol);
    Yhat:=RegCoeff(1,1);
    FOR j FROM 1 TO Order DO
      Yhat:=Yhat+RegCoeff(j+1,1)*x^j;
    END;
    Sum1:=Sum1+(Yhat-YMean)^2;
    Sum2:=Sum2+(y-YMean)^2;
  END;
  // calculate coefficient of determination
  Rsqr:=Sum1/Sum2;
  RETURN {Rsqr,RegCoeff};
END;

PowerFit

Code:
EXPORT PowerFit(MatX,VectY)
BEGIN
  LOCAL i,j;
  LOCAL lstDimX,NumRowsX,NumColsX;
  LOCAL TMatX,TVectY,RegCoeff,Rsqr;
  LOCAL YMean,Sum1,Sum2,Yhat;
  // get the number of rows and columns of matrix MatX
  lstDimX:=SIZE(MatX);
  NumRowsX:=lstDimX(1);
  NumColsX:=lstDimX(2);
  // create matrices to store transformed data
  TMatX:=MAKEMAT(1,NumRowsX,NumColsX+1);
  TVectY:=MAKEMAT(1,NumRowsX,1);
  FOR i FROM 1 TO NumRowsX DO
    TVectY(i,1):=LN(VectY(i,1));
    FOR j FROM 1 TO NumColsX DO
      TMatX(i,j+1):=LN(MatX(i,j));
    END;
  END;
  // calculate the regression coefficients
  RegCoeff:=LSQ(TMatX,TVectY);
  // calculate ymean
  Sum1:=0;
  FOR i FROM 1 TO NumRowsX DO
    Sum1:=Sum1+TVectY(i,1);
  END;
  YMean:=Sum1/NumRowsX;
  // calculate the coefficient of determination
  Sum1:=0;
  Sum2:=0;
  Yhat:=TMatX*RegCoeff;
  FOR i FROM 1 TO NumRowsX DO
    Sum1:=Sum1+(Yhat(i,1)-YMean)^2;
    Sum2:=Sum2+(TVectY(i,1)-YMean)^2;
  END;
  Rsqr:=Sum1/Sum2;
  // return the results
  RETURN {Rsqr,RegCoeff};
END;

LinearizedReg

Code:
EXPORT LinearizedReg(DSMat,lstSelXData,lstSelYData)
BEGIN
  LOCAL i,x,y;
  LOCAL lstDimX,NumRowsX;
  LOCAL MatX,VectY,RegCoeff;
  LOCAL SelXCol,ShiftX,ScaleX,PowerX;
  LOCAL SelYCol,ShiftY,ScaleY,PowerY;
  LOCAL Sum1,Sum2,YMean,Yhat,Rsqr;
  
  lstDimX:=SIZE(DSMat);
  NumRowsX:=lstDimX(1);
  // get parameters for x
  SelXCol:=lstSelXData(1);
  ShiftX:=lstSelXData(2);
  ScaleX:=lstSelXData(3);
  PowerX:=lstSelXData(4);
  // get parameters for y
  SelYCol:=lstSelYData(1);
  ShiftY:=lstSelYData(2);
  ScaleY:=lstSelYData(3);
  PowerY:=lstSelYData(4);
  // create regression matrix and vector
  MatX:=MAKEMAT(1,NumRowsX,2);
  VectY:=MAKEMAT(1,NumRowsX,1);
  // populate matrix and vector
  FOR i FROM 1 TO NumRowsX DO
    x:=ScaleX*DSMat(i,SelXCol)+ShiftX;
    y:=ScaleY*DSMat(i,SelYCol)+ShiftY;
    // transform x
    IF PowerX==0 THEN
      x:=LN(x);
    ELSE
      x:=x^PowerX;
    END;
    // transform y
    IF PowerY==0 THEN
      y:=LN(y);
    ELSE
      y:=y^PowerY;
    END;
    MatX(i,2):=x;
    VectY(i,1):=y;
  END;
  // calculate coefficient of determination
  RegCoeff:=LSQ(MatX,VectY);
  // calculate ymean
  Sum1:=0;
  FOR i FROM 1 TO NumRowsX DO
    y:=VectY(i,1);
    Sum1:=Sum1+y;
  END;
  YMean:=Sum1/NumRowsX;
  // calculate coefficient of determination
  Sum1:=0;
  Sum2:=0;
  Yhat:=MatX*RegCoeff;
  FOR i FROM 1 TO NumRowsX DO
    Sum1:=Sum1+(Yhat(i,1)-YMean)^2;
    Sum2:=Sum2+(VectY(i,1)-YMean)^2;
  END;
  Rsqr:=Sum1/Sum2;
  RETURN {Rsqr,RegCoeff};
END;

LinearizedReg3

Code:
EXPORT LinearizedReg3(DSMat,lstSelXData,lstSelZData,lstSelYData)
BEGIN
  LOCAL i,x,y,z;
  LOCAL lstDimX,NumRowsX;
  LOCAL MatX,VectY,RegCoeff;
  LOCAL SelXCol,ShiftX,ScaleX,PowerX;
  LOCAL SelYCol,ShiftY,ScaleY,PowerY;
  LOCAL SelZCol,ShiftZ,ScaleZ,PowerZ;
  LOCAL Sum1,Sum2,YMean,Yhat,Rsqr;

  lstDimX:=SIZE(DSMat);
  NumRowsX:=lstDimX(1);
  SelXCol:=lstSelXData(1);
  ShiftX:=lstSelXData(2);
  ScaleX:=lstSelXData(3);
  PowerX:=lstSelXData(4);
  // get parameters for z
  SelZCol:=lstSelZData(1);
  ShiftZ:=lstSelZData(2);
  ScaleZ:=lstSelZData(3);
  PowerZ:=lstSelZData(4);
  // get parameters for y
  SelYCol:=lstSelYData(1);
  ShiftY:=lstSelYData(2);
  ScaleY:=lstSelYData(3);
  PowerY:=lstSelYData(4);
  // create regression matrix and vector
  MatX:=MAKEMAT(1,NumRowsX,3);
  VectY:=MAKEMAT(1,NumRowsX,1);
  // populate matrix and vector
  FOR i FROM 1 TO NumRowsX DO
    x:=ScaleX*DSMat(i,SelXCol)+ShiftX;
    y:=ScaleY*DSMat(i,SelYCol)+ShiftY;
    z:=ScaleZ*DSMat(i,SelZCol)+ShiftZ;
    // transform x
    IF PowerX==0 THEN
      x:=LN(x);
    ELSE
      x:=x^PowerX;
    END;
    // transform z
    IF PowerZ==0 THEN
      z:=LN(z);
    ELSE
      z:=z^PowerZ;
    END;
    // transform y
    IF PowerY==0 THEN
      y:=LN(y);
    ELSE
      y:=y^PowerY;
    END;
    MatX(i,2):=x;
    MatX(i,3):=z;
    VectY(i,1):=y;
  END;
  // calculate regression coefficients
  RegCoeff:=LSQ(MatX,VectY);
  // calculate ymean
  Sum1:=0;
  FOR i FROM 1 TO NumRowsX DO
    y:=VectY(i,1);
    Sum1:=Sum1+y;
  END;
  YMean:=Sum1/NumRowsX;
  // calculate coefficient of determination
  Sum1:=0;
  Sum2:=0;
  Yhat:=MatX*RegCoeff;
  FOR i FROM 1 TO NumRowsX DO
    Sum1:=Sum1+(Yhat(i,1)-YMean)^2;
    Sum2:=Sum2+(VectY(i,1)-YMean)^2;
  END;
  Rsqr:=Sum1/Sum2;
  RETURN {Rsqr,RegCoeff};
END;

MLRX

Code:
EXPORT MLRX(DSMat,MaxTerms,TrnfMat)
BEGIN
  LOCAL i,j,x,Rsqr;
  LOCAL lstDimX,lstDimT,NumRowsX,NumColsX;
  LOCAL MatX,VectY,RegCoeff;
  LOCAL SelXCol,Shift,Scale,PowerX;
  LOCAL InsMat,InsCol,NumTrnf;
  LOCAL YMean,Sum1,Sum2,Yhat;

  // calculate the number of data points
  lstDimX:=SIZE(DSMat);
  NumRowsX:=lstDimX(1);
  NumColsX:=lstDimX(2);
  // get the number of transformations
  lstDimT:=SIZE(TrnfMat);
  NumTrnf:=lstDimT(1);
  // create the data matrices
  MatX:=MAKEMAT(1,NumRowsX,MaxTerms+1);
  VectY:=MAKEMAT(1,NumRowsX,1);
  FOR j FROM 1 TO NumTrnf DO
    // get the transformation/insertion parameters
    SelXCol:=TrnfMat[j,1];
    Scale:=TrnfMat[j,2];
    Shift:=TrnfMat[j,3];
    PowerX:=TrnfMat[j,4];
    InsMat:=TrnfMat[j,5];
    InsCol:=TrnfMat[j,6];
    // process all rows for current variable selection,
    // transformation, and insertion
    FOR i FROM 1 TO NumRowsX DO
      // get x
      x:=DSMat[i,SelXCol];
      // transform x by scaling and shifting
      x:=Scale*x+Shift;
      // raise x to power or take ln() value
      IF PowerX==0 THEN
        x:=LN(x);
      ELSE
        x:=x^PowerX;
      END;
      // insert in targeted matrix
      IF InsMat>0 THEN
        // insert in matrix of independent variables
        IF InsCol>(MaxTerms+1) THEN
          // display an error message
          MSGBOX("Column "+InsCol+" is outside the range of columns");
          RETURN "ERROR";
        END;
        MatX[i,InsCol]:=MatX[i,InsCol]*x;
      ELSE
        // insert in vector of dependent variable
        VectY[i,1]:=VectY[i,1]*x;
      END;
    END;
  END;
  // calculate the regression coefficients
  RegCoeff:=LSQ(MatX,VectY);
  // calculate ymean
  Sum1:=0;
  FOR i FROM 1 TO NumRowsX DO
    Sum1:=Sum1+VectY(i,1);
  END;
  YMean:=Sum1/NumRowsX;
  // calculate the correlation coefficient
  Sum1:=0;
  Sum2:=0;
  Yhat:=MatX*RegCoeff;
  FOR i FROM 1 TO NumRowsX DO
    Sum1:=Sum1+(Yhat(i,1)-YMean)^2;
    Sum2:=Sum2+(VectY(i,1)-YMean)^2;
  END;
  Rsqr:=Sum1/Sum2;
  // return the results
  RETURN {Rsqr,RegCoeff};
END;

BestLR

Code:
EXPORT BestLR(Data,lstSel,lstX,lstY)
BEGIN
  LOCAL Tx,Ty,i,j,k;
  LOCAL SelXCol,SelYCol;
  LOCAL PowerX,PowerY;
  LOCAL MatX,VectY,RegCoeff,MatRes;
  LOCAL lstDim,NumRows,ResUpdated;
  LOCAL Sum1,Sum2,Rsqr,YMean,Yhat;
  LOCAL MaxRes,NumColRes;

  // use the default transformation list
  // for variable x if lstX is empty
  IF SIZE(lstX)==0 THEN
    lstX:={-3,-2,-1,-0.5,0,0.5,1,2,3};
  END;
  // use the default transformation list
  // for variable y if lstY is empty
  IF SIZE(lstY)==0 THEN
    lstY:={-3,-2,-1,-0.5,0,0.5,1,2,3};
  END;
  // set the maximum number of results
  MaxRes:=20;
  // set the number of columns in
  // the results matrix
  NumColRes:=5;
  // get the indices for variable x and y
  SelXCol:=lstSel(1);
  SelYCol:=lstSel(2);
  // get the number of rows of data
  lstDim:=SIZE(Data);
  NumRows:=lstDim(1);
  // create the result and regression matrices
  MatX:=MAKEMAT(1,NumRows,2);
  VectY:=MAKEMAT(1,NumRows,1);
  MatRes:=    MAKEMAT(0,MaxRes,NumColRes);
  // iterate for each transformation in x
  FOR Tx FROM 1 TO SIZE(lstX) DO
    // get the current power of x
    PowerX:=lstX(Tx);
    // transform x
    IF PowerX==0 THEN
      FOR i FROM 1 TO NumRows DO
        MatX(i,2):=LN(Data(i,SelXCol));
      END;
    ELSE
      FOR i FROM 1 TO NumRows DO
        MatX(i,2):=Data(i,SelXCol)^PowerX;
      END;
    END;
    // iterate for each transformation in y
    FOR Ty FROM 1 TO SIZE(lstY) DO
      // get the current power of y
      PowerY:=lstY(Ty);
      // transform y
      IF PowerY==0 THEN
        FOR i FROM 1 TO NumRows DO
          VectY(i,1):=LN(Data(i,SelYCol));
        END;
      ELSE
        FOR i FROM 1 TO NumRows DO
          VectY(i,1):=Data(i,SelYCol)^PowerY;
        END;
      END;
      // calculate regression coefficients
      RegCoeff:=LSQ(MatX,VectY);
      // calculate ymean
      Sum1:= 0;
      FOR i FROM 1 TO NumRows DO
        Sum1:=Sum1+VectY(i,1);
      END;
      YMean:=Sum1/NumRows;
      // calculate coefficient of determination
      Sum1:=0;
      Sum2:=0;
      Yhat:=MatX*RegCoeff;
      FOR i FROM 1 TO NumRows DO
        Sum1:=Sum1+(Yhat(i,1)-YMean)^2;
        Sum2:=Sum2+(VectY(i,1)-YMean)^2;
      END;
      Rsqr:=Sum1/Sum2;
      // Rsqr is better than last entry
      // in the results matrix?
      IF Rsqr>MatRes(MaxRes,1) THEN
        ResUpdated:=0;
        // loop for each row in the
        // results matrix
        FOR i FROM 1 TO MaxRes DO
          // compare with other Rsqr values
          IF ResUpdated==0 AND Rsqr>MatRes(i,1) THEN
            // found better Rsqr for row i
            // push rows below?
            IF i<MaxRes THEN
              j:=MaxRes-1;
              REPEAT
              // FOR j FROM MaxRes-1 TO n STEP -1 DO
                FOR k FROM 1 TO NumColRes DO
                  MatRes(j+1,k):=MatRes(j,k);
                END;
                j:=j-1;
                //END;
              UNTIL j<i;
            END;
            // store new best results
            MatRes(i,1):=Rsqr;
            MatRes(i,2):=PowerY;
            MatRes(i,3):=PowerX;
            MatRes(i,4):=RegCoeff(1,1);
            MatRes(i,5):=RegCoeff(2,1);
            ResUpdated:=1;
          END;
        END;
      END;
    END; // FOR Ty
  END; // FOR Tx
  RETURN MatRes;
END;

BestMLR2

Code:
EXPORT BestMLR2(Data,lstSel,lstX,lstZ,lstY)
BEGIN
  LOCAL Tx,Tz,Ty,i,j,k;
  LOCAL SelXCol,SelZCol,SelYCol;
  LOCAL PowerX,PowerZ,PowerY;
  LOCAL MatX,VectY,RegCoeff,MatRes;
  LOCAL lstDim,NumRows,ResUpdated;
  LOCAL Sum1,Sum2,Rsqr,YMean,Yhat;
  LOCAL MaxRes,NumColRes;
  
  // use default transformation list
  // if lstX is an empty list
  IF SIZE(lstX)==0 THEN
    lstX:={-3,-2,-1,-0.5,0,0.5,1,2,3};
  END;
  // use default transformation list
  // if lstZ is an empty list
  IF SIZE(lstZ)==0 THEN
    lstZ:={-3,-2,-1,-0.5,0,0.5,1,2,3};
  END;
  // use default transformation list
  // if lstY is an empty list
  IF SIZE(lstY)==0 THEN
    lstY:={-3,-2,-1,-0.5,0,0.5,1,2,3};
  END;
  // Set the number of rows and columns in the results matrix
  MaxRes:=20;
  NumColRes:=7;
  // get the variable selectors
  SelXCol:=lstSel(1);
  SelZCol:=lstSel(2);
  SelYCol:=lstSel(3);
  // get the number of rows in matrix Data
  lstDim:=SIZE(Data);
  NumRows:=lstDim(1);
  // create the regression and results matrices
  MatX:=MAKEMAT(1,NumRows,3);
  VectY:=MAKEMAT(1,NumRows,1);
  MatRes:=MAKEMAT(0,MaxRes,NumColRes);
  // start the calculations
  // process transformations for variable x
  FOR Tx FROM 1 TO SIZE(lstX) DO
    PowerX:=lstX(Tx);
    // transform x
    IF PowerX==0 THEN
      FOR i FROM 1 TO NumRows DO
        MatX(i,2):=LN(Data(i,SelXCol));
      END;
    ELSE
      FOR i FROM 1 TO NumRows DO
        MatX(i,2):=Data(i,SelXCol)^PowerX;
      END;
    END;
    // process transformations for variable z
    FOR Tz FROM 1 TO SIZE(lstZ) DO
      PowerZ:=lstZ(Tz);
      // transform z
      IF PowerZ==0 THEN
        FOR i FROM 1 TO NumRows DO
          MatX(i,3):=LN(Data(i,SelZCol));
        END;
      ELSE
        FOR i FROM 1 TO NumRows DO
          MatX(i,3):=Data(i,SelZCol)^PowerZ;
        END;
      END;
      // process transformations for variable y
      FOR Ty FROM 1 TO SIZE(lstY) DO
        PowerY:=lstY(Ty);
        // transform y
        IF PowerY==0 THEN
          FOR i FROM 1 TO NumRows DO
            VectY(i,1):=LN(Data(i,SelYCol));
          END;
        ELSE
          FOR i FROM 1 TO NumRows DO
            VectY(i,1):=Data(i,SelYCol)^PowerY;
          END;
        END;
        // calculate regression coefficients
        RegCoeff:=LSQ(MatX,VectY);
        // calculate ymean
        Sum1:=0;
        FOR i FROM 1 TO NumRows DO
          Sum1:=Sum1+VectY(i,1);
        END;
        Ymean:=Sum1/NumRows;
        // calculate the coefficient of determination
        Sum1:=0;
        Sum2:=0;
        Yhat:=MatX*RegCoeff;
        FOR i FROM 1 TO NumRows DO
          Sum1:=Sum1+(Yhat(i,1)-Ymean)^2;
          Sum2:=Sum2+(VectY(i,1)-Ymean)^2;
        END;
        Rsqr:=Sum1/Sum2;
        // Rsqr is better than last entry
        // in the results matrix?
        IF Rsqr>MatRes(MaxRes,1) THEN
          ResUpdated:=0;
          // check which row to insert better
          // regression results
          FOR i FROM 1 TO MaxRes DO
            // insert new results in row i?
            IF ResUpdated==0 AND Rsqr>MatRes(i,1) THEN
              // inserting inside the results matrix?
              IF i<MaxRes THEN
                // downward copy rows from MaxRes-1 to i
                j:=MaxRes-1;
                REPEAT
                  FOR k FROM 1 TO NumColRes DO
                    MatRes(j+1,k):=MatRes(j,k);
                  END;
                  j:=j-1;
                UNTIL j<i;
              END;
              // insert better results in row i
              MatRes(i,1):=Rsqr;
              MatRes(i,2):=PowerY;
              MatRes(i,3):=PowerX;
              MatRes(i,4):=PowerZ;
              FOR k FROM 1 TO 3 DO
                MatRes(i,4+k):=RegCoeff(k,1);
              END;
              ResUpdated:=1;
            END;
          END;
        END;
      END; // FOR Ty
    END; // FOR Tz
  END; // FOR Tx
  RETURN MatRes;
END;

BestPolyReg

Code:
EXPORT BestPolyReg(Data,SelXCol,SelYCol,MinOrder,MaxOrder)
BEGIN
  LOCAL i,k,x,y,Order;
  LOCAL MatX,VectY;
  LOCAL lstDim,NumRows;
  LOCAL Sum1,Sum2,Sum3,AICc,BestAICc;
  LOCAL Rsqr,RsqrAdj,YMean,Yhat,RegCoeff;
  LOCAL BestOrder,BestRsqr,BestRsqrAdj,BestRegCoeff;
  lstDim:=SIZE(Data);
  NumRows:=lstDim(1);
  IF MinOrder<1 THEN
    MinOrder:=1;
  END;
  // initialize best regression data
  BestOrder=0;
  BestRsqr:=0;
  BestRsqrAdj=0;
  BestAICc:=1E499;
  // initialize vector y
  VectY:=MAKEMAT(1,NumRows,1);
  // calculate ymean … needs to be done once!
  Sum1:=0;
  FOR i FROM 1 TO NumRows DO
    y:=Data(i,SelYCol);
    VectY(i,1):=y;
    Sum1:=Sum1+y;
  END;
  YMean:=Sum1/NumRows;
  // calculate sum of y – ymean squared
  Sum2:=0;
  FOR i FROM 1 TO NumRows DO
    y:=Data(i,SelYCol);
    Sum2:=Sum2+(y-YMean)^2;
  END;
  // iterate for the specified range of polynomial orders
  FOR Order FROM MinOrder TO MaxOrder DO
    // (re)create the matrix MatX
    MatX:=MAKEMAT(1,NumRows,1+Order);
    // fill the columns to to Order+1 with x^power values
    FOR i FROM 1 TO NumRows DO
      x:=Data(i,SelXCol);
      FOR k FROM 1 TO Order DO
        MatX(i,k+1):=x^k;
      END;
    END;
    // calculate regression coefficients
    RegCoeff:=LSQ(MatX,VectY);
    // calculate the coefficient of determination
    Sum1:=0;
    Sum3:=0;
    Yhat:=MatX*RegCoeff;
    FOR i FROM 1 TO NumRows DO
      Sum1:=Sum1+(Yhat(i,1)-YMean)^2;
      Sum3:=Sum3+(VectY(i,1)-Yhat(i,1))^2;
    END;
    Rsqr:=Sum1/Sum2;
    // calculate the adjusted coefficient of determination
    RsqrAdj:=1-(1-Rsqr)*(NumRows-1)/(NumRows-Order-1);
    k:=Order+1;
    // if Sum3 is 0 then adjust it to a small value
    // to avoid getting a LN(0) error
    IF Sum3==0 THEN
      Sum3:=1E-499;
    END;
    // calculate AICc statistic
    AICc:=NumRows*LN(Sum3/NumRows)+2*k+(2*k*(k+1))/(NumRows-k-1);
    // display intermediate results
    PRINT("Order="+Order+", AICc="+AICc+", Sum3="+Sum3);
    // found a better fit?
    IF AICc<BestAICc THEN
      // update best regression data
      BestOrder:=Order;
      BestRsqr:=Rsqr;
      BestRsqrAdj:=RsqrAdj;
      BestAICc:=AICc;
      BestRegCoeff:=RegCoeff;
    END;
  END;
  RETURN {BestOrder,BestRsqr,BestRsqrAdj,BestAICc,BestRegCoeff};
END;

BestMLR3

Code:
EXPORT BestMLR3(Data,lstSel,lstX,lstZ,lstT,lstY)
BEGIN
  LOCAL Tx,Tz,Ty,Tt,i,j,k;
  LOCAL SelXCol,SelZCol,SelTCol,SelYCol;
  LOCAL PowerX,PowerZ,PowerT,PowerY;
  LOCAL MatX,VectY,RegCoeff,MatRes;
  LOCAL lstDim,NumRows,ResUpdated;
  LOCAL Sum1,Sum2,Rsqr,YMean,Yhat;
  LOCAL MaxRes,NumColRes;

  // use default transformation list
  // if lstX is an empty list
  IF SIZE(lstX)==0 THEN
    lstX:={-3,-2,-1,-0.5,0,0.5,1,2,3};
  END;
  // use default transformation list
  // if lstZ is an empty list
  IF SIZE(lstZ)==0 THEN
    lstZ:={-3,-2,-1,-0.5,0,0.5,1,2,3};
  END;
  // use default transformation list
  // if lstT is an empty list
  IF SIZE(lstT)==0 THEN
    lstT:={-3,-2,-1,-0.5,0,0.5,1,2,3};
  END;
  // use default transformation list
  // if lstY is an empty list
  IF SIZE(lstY)==0 THEN
    lstY:={-3,-2,-1,-0.5,0,0.5,1,2,3};
  END;
  MaxRes:=20;
  NumColRes:=9;
  // get the variable selectors
  SelXCol:=lstSel(1);
  SelZCol:=lstSel(2);
  SelTCol:=lstSel(3);
  SelYCol:=lstSel(4);
  // get the number of rows in matrix Data
  lstDim:=SIZE(Data);
  NumRows:=lstDim(1);
  // create the regression and results matrices
  MatX:=MAKEMAT(1,NumRows,3);
  VectY:=MAKEMAT(1,NumRows,1);
  MatRes:=MAKEMAT(0,MaxRes,NumColRes);
  // start the calculations
  // process transformations for variable x
  FOR Tx FROM 1 TO SIZE(lstX) DO
    PowerX:=lstX(Tx);
    // transform x
    IF PowerX==0 THEN
      FOR i FROM 1 TO NumRows DO
        MatX(i,2):=LN(Data(i,SelXCol));
      END;
    ELSE
      FOR i FROM 1 TO NumRows DO
        MatX(i,2):=Data(i,SelXCol)^PowerX;
      END;
    END;
    // process transformations for variable z
    FOR Tz FROM 1 TO SIZE(lstZ) DO
      PowerZ:=lstZ(Tz);
      // transform z
      IF PowerZ==0 THEN
        FOR i FROM 1 TO NumRows DO
          MatX(i,3):=LN(Data(i,SelZCol));
        END;
      ELSE
        FOR i FROM 1 TO NumRows DO
          MatX(i,3):=Data(i,SelZCol)^PowerZ;
        END;
      END;
      // process transformations for variable t
      FOR Tt FROM 1 TO SIZE(lstT) DO
        PowerT:=lstT(Tt);
        // transform t
        IF PowerT==0 THEN
          FOR i FROM 1 TO NumRows DO
            MatX(i,4):=LN(Data(i,SelTCol));
          END;
        ELSE
          FOR i FROM 1 TO NumRows DO
            MatX(i,4):=Data(i,SelTCol)^PowerT;
          END;
        END;
        // process transformations for variable y
        FOR Ty FROM 1 TO SIZE(lstY) DO
          PowerY:=lstY(Ty);
          // transform y
          IF PowerY==0 THEN
            FOR i FROM 1 TO NumRows DO
              VectY(i,1):=LN(Data(i,SelYCol));
            END;
          ELSE
            FOR i FROM 1 TO NumRows DO
              VectY(i,1):=Data(i,SelYCol)^PowerY;
            END;
          END;
          // calculate regression coefficients
          RegCoeff:=LSQ(MatX,VectY);
          // calculate ymean
          Sum1:=0;
          FOR i FROM 1 TO NumRows DO
            Sum1:=Sum1+VectY(i,1);
          END;
          YMean:=Sum1/NumRows;
          // calculate the coefficient of determination
          Sum1:=0;
          Sum2:=0;
          Yhat:=MatX*RegCoeff;
          FOR i FROM 1 TO NumRows DO
            Sum1:=Sum1+(Yhat(i,1)-YMean)^2;
            Sum2:=Sum2+(VectY(i,1)-YMean)^2;
          END;
          Rsqr:=Sum1/Sum2;
          // Rsqr is better than last entry
          // in the results matrix?
          IF Rsqr>MatRes(MaxRes,1) THEN
            ResUpdated:=0;
            // check which row to insert better
            // regression results
            FOR i FROM 1 TO MaxRes DO
              // insert new results in row i?
              IF ResUpdated==0 AND Rsqr>MatRes(i,1) THEN
                // inserting inside the results matrix?
                IF i<MaxRes THEN
                  // downward copy rows from MaxRes-1 to i
                  j:=MaxRes-1;
                  REPEAT
                    FOR k FROM 1 TO NumColRes DO
                      MatRes(j+1,k):=MatRes(j,k);
                    END;
                    j:=j-1;
                  UNTIL j<i;
                END;
                // insert better results in row i
                MatRes(i,1):=Rsqr;
                MatRes(i,2):=PowerY;
                MatRes(i,3):=PowerX;
                MatRes(i,4):=PowerZ;
                MatRes(i,5):=PowerT;
                FOR k FROM 1 TO 4 DO
                  MatRes(i,5+k):=RegCoeff(k,1);
                END;
                ResUpdated:=1;
              END;
            END;
          END;
        END; // FOR Ty
      END; // FOR Tt
    END; // FOR Tz
  END; // FOR Tx
  RETURN MatRes;
END;

RErMLR2

Code:
EXPORT RErMLR2(MatX,VectY)
BEGIN
  LOCAL i,j,y,Rsqr;
  LOCAL lstDimX,NumRows,NumCols;
  LOCAL YMean,Sum1,Sum2,Yhat;
  LOCAL TMatX,VectOnes,RegCoeff;
  // calculate the number of observations and variables
  lstDimX:=SIZE(MatX);
  NumRows:=lstDimX(1);
  NumCols:=lstDimX(2);
  // create transformation matrices
  TMatX:=MAKEMAT(1,NumRows,NumCols+1);
  VectOnes:=MAKEMAT(1,NumRows,1);
  FOR i FROM 1 TO NumRows DO
    y:=VectY(i,1);
    TMatX(i,1):=1/y;
    FOR j FROM 1 TO NumCols DO
      TMatX(i,j+1):=MatX(i,j)/y;
    END;
  END;
  // calculate the regression coefficients
  RegCoeff:=LSQ(TMatX,VectOnes);
  // calculate ymean
  Sum1:=0;
  FOR i FROM 1 TO NumRows DO
    y:=VectY(i,1);
    Sum1:=Sum1+y;
  END;
  YMean:=Sum1/NumRows;
  // calculate the coefficient of determination
  Sum1:=0;
  Sum2:=0;
  FOR i FROM 1 TO NumRows DO
    y:=VectY(i,1);
    Yhat:=RegCoeff(1,1);
    FOR j FROM 1 TO NumCols DO
      Yhat:=Yhat+RegCoeff(j+1,1)*MatX(i,j);
    END;
    Sum1:=Sum1+((Yhat-YMean)/y)^2;
    Sum2:=Sum2+((y-YMean)/y)^2;
  END;
  Rsqr:=Sum1/Sum2;
  // return the results as a list
  RETURN {Rsqr,RegCoeff};
END;

RErPolyReg

Code:
EXPORT RErPolyReg(DSMat,SelXCol,SelYCol,Order)
BEGIN
  LOCAL i,j,x,y;
  LOCAL lstDimX,NumRows;
  LOCAL MatX,VectOnes,RegCoeff;
  LOCAL YMean,Sum1,Sum2,Yhat,Rsqr;

  // get the number of rows
  lstDimX:=SIZE(DSMat);
  NumRows:=lstDimX(1);
  // create the regression matrices
  MatX:=MAKEMAT(1,NumRows,Order+1);
  VectOnes:=MAKEMAT(1,NumRows,1);
  // populate the regression matrices
  FOR i FROM 1 TO NumRows DO
    y:=DSMat(i,SelYCol);
    x:=DSMat(i,SelXCol);
    MatX(i,1):=1/y;
    FOR j FROM 1 TO Order DO
      MatX(i,j+1):=(x^j)/y;
    END;
  END;
  // calculate the regression coefficients
  RegCoeff:=LSQ(MatX,VectOnes);
  // calculate ymean
  Sum1:=0;
  FOR i FROM 1 TO NumRows DO
    Sum1:=Sum1+DSMat(i,SelYCol);
  END;
  YMean:=Sum1/NumRows;
  // calculate the coefficient of determination
  Sum1:=0;
  Sum2:=0;
  FOR i FROM 1 TO NumRows DO
    x:=DSMat(i,SelXCol);
    y:=DSMat(i,SelYCol);
    Yhat:=RegCoeff(1,1);
    FOR j FROM 1 TO Order DO
      Yhat:=Yhat+RegCoeff(j+1,1)*x^j;
    END;
    Sum1:=Sum1+((Yhat-YMean)/y)^2;
    Sum2:=Sum2+((y-YMean)/y)^2;
  END;
  Rsqr:=Sum1/Sum2;
  // return the list of results
  RETURN {Rsqr,RegCoeff};
END;

RErPowerFit

Code:
EXPORT RErPowerFit(MatX,VectY)
BEGIN
  LOCAL i,j,y,LnY;
  LOCAL lstDimX,NumRows,NumCols;
  LOCAL TMatX,VectOnes,RegCoeff,Rsqr;
  LOCAL YMean,Sum1,Sum2,Yhat;

  lstDimX:=SIZE(MatX);
  NumRows:=lstDimX(1);
  NumCols:=lstDimX(2);
  TMatX:=MAKEMAT(1,NumRows,NumCols+1);
  VectOnes:=MAKEMAT(1,NumRows,1);
  Sum1:=0;
  FOR i FROM 1 TO NumRows DO
    LnY:=LN(VectY(i,1));
    Sum1:=Sum1+LnY;
    TMatX(i,1):=1/LnY;
    FOR j FROM 1 TO NumCols DO
      TMatX(i,j+1):=LN(MatX(i,j))/LnY;
    END;
  END;
  // calculate regression coefficients
  RegCoeff:=LSQ(TMatX,VectOnes);
  // calculate mean ln(y)
  YMean:=Sum1/NumRows;
  // calculate coefficient of determination
  Sum1:=0;
  Sum2:=0;
  FOR i FROM 1 TO NumRows DO
    Yhat:=RegCoeff(1,1);
    FOR j FROM 1 TO NumCols DO
      Yhat:=Yhat+RegCoeff(j+1,1)*LN(MatX(i,j));
    END;
    y:=LN(VectY(i,1));
    Sum1:=Sum1+((Yhat-YMean)/y)^2;
    Sum2:=Sum2+((y-YMean)/y)^2;
  END;
  Rsqr:=Sum1/Sum2;
  RETURN {Rsqr,RegCoeff};
END;
Good info, thanks for sharing!

A quick and dirty linear regression can always be set up with LSQ, but the programs are helpful as they also contain regression diagnostics.

Is there a repository for other HP 39gII programs?
(09-04-2014 03:57 PM)Helge Gabert Wrote: [ -> ]Is there a repository for other HP 39gII programs?

Not specifically for the HP 39gII. It would be in the general area. I will copy my reply to the general area.
my first draft for ANOVA with Prime is here.
(09-04-2014 03:28 PM)Namir Wrote: [ -> ]BestMLR2

Code:
EXPORT BestMLR2(Data,lstSel,lstX,lstZ,lstY)
BEGIN
  ...
  LOCAL Sum1,Sum2,Rsqr,YMean,Yhat;
  ...
        // calculate ymean
        Sum1:=0;
        FOR i FROM 1 TO NumRows DO
          Sum1:=Sum1+VectY(i,1);
        END;
        Ymean:=Sum1/NumRows;
        // calculate the coefficient of determination
        Sum1:=0;
        Sum2:=0;
        Yhat:=MatX*RegCoeff;
        FOR i FROM 1 TO NumRows DO
          Sum1:=Sum1+(Yhat(i,1)-Ymean)^2;
          Sum2:=Sum2+(VectY(i,1)-Ymean)^2;
...
END;


hi Namir,
in BestMLR2 there is a little typo: YMean and Ymean; to correct it I simply change the namevar in Local definition...
However is a very useful program!
All the suite of these programs should be inside the FW of Prime, perhaps with other programs to analyze multiple regression, thank you to have done them!

Salvo
Reference URL's