Coeficienții ecuației de regresie simplă

Modelul de regresie simplă

pentru calculul coeficienților a și b se folosesc formulele 

subprogramele de calcul pentru coeficienții a și b sunt următoarele:

SUMVEC - sumă elemente ale unui vector S=X(1)+X(2)+....+X(N);

PATRAT - ridicare la pătrat a elementelor unui vector Z(I)=X(I)*X(I);

PRODEL - calcul produse elemente ale unui vector Z(I)=X(I)*Y(I);

COEF - calculează coeficineții A și B după formulele de mai sus.

         SUBROUTINE SUMVEC(X,N,S)

         S=0.

         DO   10  I=1,N

         S=S+X(I)

10     CONTINUE

         RETURN

         END

         SUBROUTINE PATRAT(X,N,S2)

         S2=0.

         DO   10  I=1,N

         S2=S2+X(I)*X(I)

10     CONTINUE

         RETURN

         END

         SUBROUTINE PRODEL(X,Y,Z,N)

         DO   10  I=1,N

         Z(I)=X(I)*Y(I)

10     CONTINUE

         RETURN

         END

         SUBROUTINE COEF(X,Y,N,A,B,IK)

         IK=0

         IF(N.LE.0) GO TO 10

         AN=N

         IK=1

         CALL PRODEL(X,Y,XY,N)

         CALL  SUMVEC(X,N,SX)

         CALL  SUMVEC(XY,N,SXY)

         CALL  SUMVEC(Y,N,SY)

         CALL  PATRAT(X,XN,SXX)

         A=(AN*SXY-SX*SY)/(AN*SXX-SX*SX)

         B=(SY-A*SX)/AN

10     CONTINUE

         RETURN

         END

Variabila IK este 0 dacă nu se calculează coeficienții că seriile de date sunt vide și IK=1 dacă acești coeficienți se calculează. Subprogramul COEF() se optimizează și dacă se elimină puzderia de apeluri de subprograme, dar devine mai mărișor și cu instrucțiuni ceva mai regrupate după cum urmează:

         SUBROUTINE COEF(X,Y,N,A,B,IK)

         IK=0

         IF(N.LE.0) GO TO 20

         AN=N

         IK=1

         SX=0.

         SY=0.

         SXX=0.

         SXY=0.

         DO   10  I=1,N

         SX=SX+X(I)

         SY=SY+Y(I)

         SXX=SXX+X(I)*X(I)

         SXY=SXY+X(I)*Y(I)

10     CONTINUE

         A=(AN*SXY-SX*SY)/(AN*SXX-SX*SX)

         B=(SY-A*SX)/AN

20     CONTINUE

         RETURN

         END

Dacă subprogramul ar fi scris în limbajul C++ cu siguranță că inițializările s-ar fi scris SX=SY=SXX=SXY=0.; ceea ce ar fi fost cu mult mai drăguț. Coeficienții ecuației de regresie simplă se calculează fie apelând tot felul de alte subprograme simple, fie scriind un subprogram ce include toate calculele. Trebuie ales un echilibru între a apela subprograme și a avea redundanță în program prin nereutilizare de componente de bibliotecă. Programatorul e suveran în programul pe care-l scrie și el alege dacă vrea să fie genial sau penibil, riscuș fiind în ambele cazuri 50% și 50%.

(29 decembrie 2017)

Alegerea identificatorilor

În limbajul FORTRAN IV căci cu el se scriau programe până în 1977 numele de variabile care se construiau trebuiau să respecte următoarele reguuli:

- lungimea numelui de variabilă era maxim de 6 caractere;

- caracterele utilizate sunt literele mari A-Z, cifrele 0-9 și semnul $;

- primul caracter al numelui de variabilă este obligatoriu o literă.

Astfel, construcțiile:

VALOARE

COST

PROFIT

CERC

G981

$VAR

SUMA3

sunt nume de variabile corecte.

Construcțiile:

1RAZA   nu este corectă că are prim caracter altceva decât o literă;

RAZA+  nu este corectă că are simbolul + care nu este nici literă și nici cifră;

PR_PL    nu este corectă pentru că are în interior caracterul underscore cod HTML _;

A1234567B nu este corectă pentru că are mai mult de 6 caractere în componență.

Chiar dacă erau restricții așa de severe, programatorii de talent știau să aleagă astfel de denumiri ale variabilelor astfel încât să asigure un nivel de ortogonalitate pentru a nu le confunda unele cu altele. Se făceau o serie de convenții precum:

- variabila I era folosită pentru traversarea elementelor unui șir sau a liniilor matricei;

- variabila J era folosită pentru traversarea elementelor de pe coloanele unei matrice;

- variabila TEMP era folosită pentru rezultate temporare, ca și variabila AUX;

- variabilele M și N erau utilizate pentru a stabili număr de elemente sau de linii sau coloane;

- variabila XMED era folosită pentru media aritmetică;

- variabila DISP era folosită pentru dispersie;

- postfixul SORT era folosit pentru a marca faptul că rezultatul este sortat XSORT, YSORT;

- prefixul STOC era folosit pentru stocuri STOCPL, STOCRE, STOCEF;

- postfixul PL era utilizat pentru planificat, RE pentru realizat și EF pentru efectiv.

Fiecare programator își construia propriile reguli de a defini nume de variabile, dar ideia era ca programul să fie cât mai lizibil și el ca programator să nu depună eforturi inutile nici la depanare și nici la mentenanță, ci să meargă la sigur, având reguli pe care și le-a definit și le-a și respectat. 

Programatorii care au lucrat foarte ermtic, definind nume de variabile după niște reguli foarte codificate au avut mari probleme în a depana programele, căci este foarte dificil să se umble într-o secvență care pleacă de la o relație matematică de forma:

să fie recunoscută cu lejeritate în construcția:

S234=S65*S987*S987

în loc de a scrie:

SCERC=PI*R*R

Este rezonabil mai ales când se construiesc subprograme să se utilizeze modalități omogene de definire a numelor de variabile EPS pentru acel epsilon care controlează precizia, ITER pentru numărul de iterații și pentru variabila de stare IK sau ISTARE sau pentru codul de eroare IERR și toți programatorii vor folosi fără retitență bibliotecile de subprograme, altfel le vor evita că sunt greoaie. Alegerea identificatorilor intră în arta programatorilor și adevărații programatori au și reguli de a construi nume de variabile care săă facă programele livibile și ușor de manevrat, fără ca expresiile să fie groaznic de sofusticate. În FORTRAN IV nu există cuvinte cheie și programatorii ar avea la dispoziție să folosească DO, IF, READ, GO TO RETURN ca nume de variabile. M-am jucat și eu cândva cu acestă facilitate și noroc că a fost doar o joacă pentru că un program cu astfel de variabile este o calamitate la depanare, căci există riscul de a distruge din prea mult zel chiar instrucțiuni corect scrise.

(28 decembrie 2017)

Activarea și dezactivarea de termeni

În informatică nu este rezonabil să fie șterse nivi articole din fișiere și nici termeni ai seriilor de date. Ceea ce ar trebui șters, trebuie dezactivat. În acest sens se asociază fiecărui element din șir un alt element care are valoarea 0 dacă este considerat neactiv sau șters din punct de vedere logic și nu fizic. La fișiere în articol se definește un câmp de un byte care este 0 dacă articolul este neactiv sau șters logic și 1 dacă este activ. Prin testarea acelui câmp, dacă are valoarea 1 se for face prelucrări, iar sacă este zero nu va fi utilizat nici articolul și nici elementul din șir.

Subprogramul care însumează elementele active ale șirului X de N componente va avea un subprogram care mai primește ca parametru un vector LL tot de N componente și LL(I) =0 ceea  ce înseamnă că elementul nu este activ, este șters logic și nu intră în calcule, iar dacă LL(I)=1 elementul X(I) este activ și se folosește în prelucrări.

             SUBROUTINE SUMVEC (X,LL,N,SUM)

             SUM=0.

             DO  10 I=1,N

             IF(LL(I).EQ.0) GO TO 10

             SUM=SUM+X(I)

10         CONTINUE

             RETURN

             END

Activând și dezactivând termeni se construiesc tot felul de eșantioane și se verifică tot felul de ipoteze. Prin operații simple se procedează la reconstituirea de toturi de date fără a fi nevoie să se refacă nici seriile de date și nici fișiere. Acum în prelucrările moderne se lucrează în mod curent cu vectori sau câmpuri de stare fără a mai se proceda șa ștergeri de elemente sau la ștergeri de articole.

(28 decembrie 2017)

Inversa matricei nesingulare

Inversarea matricei nesingulare este subprogramul cel mai căutat căci nenumărate aplicații folosesc inversa matricei, de la ecuația de regresie până la rezolvarea de sisteme liniare.

Am tot căutat subprograme pentru inversarea unei matrice și dintre toate am găsit unul interesant în cartea PROBLEME DE ANALIZĂ NUMERICĂ REZOLVATE CU CALCULATORUL avându-i ca autori pe acad. Gheorghe MARINESCU, Irona RIZZOLI, Ileana POPESCU și Cristina ȘTEFAN, carte publicată îmn Editura Academiei RSR, în anul 1987 și care la paginile 171-173 prezintă metoda iterativă HOTELLING-BODEWING de inversare a unei matrice. Se construiește șirul lui NEWTON de torma X(n+1) = X(n) +(I-X(n)*A)*X(n), unde A, I, X, sunt matrice pătrate cu m linii și m coloane, iar cifrele cuprinse între paranteze acolo la indici arată numărul iterației. X(n+1) înseamnă matricea X la iterația n+1, iar X(n) înseamnă matricea X la iterația n. I este matricea unitate, A este matricea ce trebuie inversată. Subprogramul are ca prim pas evaluarea expresiei X1=X0 +(I+X0*A)*X0. În final matricea X1 va conține matricea inversă, iar X0 este matricea de start care se obține prin convenție dintr-o evaluare de forma: X0 = 1./(maxaij) * I. Variabila maxaij este de fapt elementul maxim dintre elementele matricei A. În subprogramul INVMAT, X0 se va referi la matricea de la iterația precedentă și X1 se va referi la matyricea de la iterația următoare.

Se observă că aici trebuie să fie apelate alte subprograme precum:

- subprogram de copiere matrice;

- subprogram de adunare matrice;

- subprogram de scădere matrice;

- subprogram de calcul normă matrice;

- subprogram de înmulțire matrice.

Față de programul dat în cartea pe care am citat-o am făcut o serie de modificări, căci m-am gândit ca utilizatorul să fie scutit în a da el acolo valoarea de start X0, eroarea, criteriu de stopare a procesului iterativ și matricea unitate.

         SUBROUTINE START(A,N,X0)

         MAXAIJ=A(1,1)

         DO  10 I=1,N

         DO  10 J=1,N

         IF(MAXAIJ.GT.A(I,J)) GO TO 10

         MAXAIJ=A(I,J)

10     CONTINUE

         DO  20 I=1,N

         DO  20 J=1,N

20     X0(I,J)=0.

         ALFA=1./MAXAIJ

         DO  30 I=1,N

30     X0(I,I)=ALFA

         RETURN

         END

Această subrutină construiește matricea de start X0 ca fiind matricea unitate înmulțită cu inversul elementului maxim din matricea A.

Pentru scăderea a două matrice pătrate având n linii și n coloane se foolosește subprogramul:

         SUBROUTINE ADMAT(A,B,C,N)

         MAXAIJ=A(1,1)

         DO  10 I=1,N

         DO  10 J=1,N

         C(I,J)=A(I,J)+B(I,J) 

10     CONTINUE

         RETURN

         END

Pentru scăderea a două matrice pătrate având n linii și n coloane se foolosește subprogramul:

         SUBROUTINE ADMAT(A,B,C,N)

         MAXAIJ=A(1,1)

         DO  10 I=1,N

         DO  10 J=1,N

         C(I,J)=A(I,J)-B(I,J) 

10     CONTINUE

         RETURN

         END

Pentru copierea matricei A în matricea C, ambele având n linii și n coloane se foolosește subprogramul:

         SUBROUTINE SCMAT(A,C,N)

         MAXAIJ=A(1,1)

         DO  10 I=1,N

         DO  10 J=1,N

         C(I,J)=A(I,J)

10     CONTINUE

         RETURN

         END

Generarea unei matrice unitate de ordin n se face cu subprogramul:

         SUBROUTINE UNIMT(N,UNIT)

         DO  20 I=1,N

         DO  10 J=1,N

         UNIT(I,J)=0.

10     CONTINUE

         UNIT(I,I)=1.

20     CONTINUE

         RETURN

         END

Norma ANOR a matricei A este calculată cu subprogramul:

         SUBROUTINE NMMAT(A,N,ANOR)

         ANOR=ABS(A1,1))

         DO  20 I=1,N

         DO  10 J=1,N

         IF(ANOR.GT.ABS(A(I,J)) GO TO 10

         ANOR=ABS(A(I,J))

10     CONTINUE

         RETURN

         END

Pentru înmulțirea a două matrice A și B pătrate având n linii și n coloane, rezultatul fiind matricea C tot cu n linii și n coloane se foolosește subprogramul:

         SUBROUTINE PRMAT(A,C,N)

         MAXAIJ=A(1,1)

         DO  20 I=1,N

         DO  20 J=1,N

         CIJ=0.

         DO  10 K=1,N

         CIJ=CIJ+A(I,K)*B(K,J)

10     CONTINUE

         C(I,J)=CIJ

20     CONTINUE

         RETURN

         END

Subprogramul pentru calculul matricei inverse AINV a matricei A cu n linii și n coloane este:

         SUBROUTINE INVMAT(A,AINV,N,IK)

         ITER=0

         Q=.01

         EPS=.0001

         CALL START(A,N,X0)

10     CALL UNIMT(N,U)

         CALL PRMAT(A,X0,X1,N)

         CALL SCMAT(UNIT,X1,DIF,N)

         ITER=ITER+1

         CALL PRMAT(DIF,X0,X1,N)

         CALL ADMAT(X0,X1,X1,N)

         CALL NMMAT(X0,N,ANOR)

         IF(Q**ITER*ANOR/(1.-Q).LE.EPS) GO TO 20

         CALL CPMAT(X1,X0,N)        

         GO TO 10

20     CONTINUE

         CALL CPMAT(X1,AINV,N)

         RETURN

         END

Față de varianta inițială, aici sunt date valorile care dau finețea rezultatului, fără a da flexibilitate utilizatorului să obțină el precizia dorită, din dorința de a nu avea mulți paramertrii în listă.

(27 decembrie 2017)