GGE SAS code - Githubissues

ArminSaed / GGE_SAS

Calculating GGE, genotype plus genotype by environment interaction, for Plant Breeding purposes. the code is written by Dr. Armin Saed-Moucheshi

0 stars 0 forks source link

/by Dr. A. Saed-Moucheshi/

DATA RAW;

INPUT ENV$ GEN$ YLD; yld=yld/10; cards; 1 1 8.1724 1 2 7.4259 1 3 5.7265 10 1 6.7013 10 2 4.5092 10 3 5.4628 11 1 4.2235 11 2 3.7342 11 3 5.0893 12 1 5.9093 12 2 6.3898 12 3 6.3797 13 1 6.948 13 2 6.1017 13 3 5.8165 14 1 8.9242 14 2 6.4404 14 3 5.9076 15 1 8.0515 15 2 5.4309 15 3 4.1781 16 1 5.8129 16 2 5.1571 16 3 4.4935 17 1 6.0399 17 2 4.8329 17 3 3.6025 18 1 8.3884 18 2 6.6273 18 3 4.9391 19 1 8.1866 19 2 4.8112 19 3 3.2301 2 1 7.5141 2 2 5.0989 2 3 5.3435 20 1 6.9851 20 2 6.255 20 3 5.0124 21 1 5.0792 21 2 2.9249 21 3 3.0342 22 1 8.2205 22 2 5.6061 22 3 5.9177 23 1 4.9317 23 2 4.9865 23 3 5.0533 24 1 6.4893 24 2 5.1367 24 3 4.2885 25 1 8.5063 25 2 5.9433 25 3 4.3384 26 1 8.1424 26 2 6.0032 26 3 4.6498 27 1 7.8328 27 2 5.4007 27 3 4.5176 28 1 6.1548 28 2 5.9983 28 3 4.3815 29 1 0.1252 29 2 6.1298 29 3 6.0535 3 1 5.8876 3 2 5.0207 3 3 4.2272 4 1 6.5973 4 2 5.5915 4 3 4.588 5 1 8.8154 5 2 7.2745 5 3 4.4955 6 1 6.6465 6 2 4.8949 6 3 4.6416 7 1 6.2611 7 2 3.59 7 3 4.0908 8 1 7.2591 8 2 5.1967 8 3 5.3745 9 1 6.7069 9 2 5.6807 9 3 4.4745

;

/ analysis linear-bilinear model / PROC GLM DATA=RAW OUTSTAT=STATS ; CLASS ENV GEN; MODEL YLD = ENV GEN ENVGEN/SS4; DATA STATS2; SET STATS ; DROP NAME TYPE; IF SOURCE = 'ERROR' THEN DELETE; / values obtained from previous analysis / MSE=2.27547331; DFE=168; NREP=3; SS=SSNREP; MS=SS/DF; F=MS/MSE; PROB=1-PROBF(F,DF,DFE); PROC PRINT DATA=STATS2 NOOBS; VAR SOURCE DF SS MS F PROB; / define AMMI model / PROC GLM DATA=RAW NOPRINT; CLASS ENV GEN; MODEL YLD = ENV GEN / SS4 ; OUTPUT OUT=OUTRES R=RESID; PROC SORT DATA=OUTRES; BY GEN ENV; PROC TRANSPOSE DATA=OUTRES OUT=OUTRES2; BY GEN; ID ENV; VAR RESID; PROC IML; USE OUTRES2; READ ALL INTO RESID; NGEN=NROW(RESID); NENV=NCOL(RESID); USE STATS2; READ VAR {MSE} INTO MSEM; READ VAR {DFE} INTO DFEM; READ VAR {NREP} INTO NREP; CALL SVD (U,L,V,RESID); MINIMO=MIN(NGEN,NENV); L=L[1:MINIMO,]; SS=(L##2)NREP; SUMA=SUM(SS); PERCENT=((1/SUMA)#SS)100; MINIMO=MIN(NGEN,NENV); PERCENTA=0; DO I = 1 TO MINIMO; DF=(NGEN-1)+(NENV-1)-(2I-1); DFA=DFA//DF; PORCEACU=PERCENT[I,]; PERCENTA=PERCENTA+PORCEACU; PERCENAC=PERCENAC//PERCENTA; END; DFE=J(MINIMO,1,DFEM); MSE=J(MINIMO,1,MSEM); SSDF=SS||PERCENT||PERCENAC||DFA||DFE||MSE; L12=L##0.5; SCOREG1=U[,1]#L12[1,]; SCOREG2=U[,2]#L12[2,]; SCOREG3=U[,3]#L12[3,]; SCOREE1=V[,1]#L12[1,]; SCOREE2=V[,2]#L12[2,]; SCOREE3=V[,3]#L12[3,]; SCOREG=SCOREG1||SCOREG2||SCOREG3; SCOREE=SCOREE1||SCOREE2||SCOREE3; SCORES=SCOREG//SCOREE; CREATE SUMAS FROM SSDF; APPEND FROM SSDF; CLOSE SUMAS; CREATE SCORES FROM SCORES; APPEND FROM SCORES ; CLOSE SCORES; / obtaining the biplot’s polygon and its perpendiculars / d1=scoreg[,1:2][cvexhull(scoreg[,1:2])[loc(cvexhull(scoreg[,1:2])>0),],]; d=d1//d1[1,]; xxx=J(nrow(d)-1,1,0); yyy=J(nrow(d)-1,1,0); ppp={0 1,1 0}; do i=1 to nrow(d)-1 ; dd=d[i:i+1,]; if dd[1,1]>dd[2,1] then ddd=pppdd; else ddd=dd; p=(ddd[2,2]-ddd[1,2])/(ddd[2,1]-ddd[1,1]) ;

if p<0 then ss=1 ; else ss=-1 ; r=tan((180-90- abs(atan(p)180/3.14156))3.14156/180)ss ; aa=(ddd[1,2]+ddd[2,2])/2-p(ddd[1,1]+ddd[2,1])/2; xx=aa/(r-p) ; if abs(r)<1 then xxx[i,]=1; else xxx[i,]=1/abs(r); if xx<0 then xxx[i,]=-xxx[i,] ; else xxx[i,]=xxx[i,]; yyy[i,]=xxx[i,]r; end; kk=xxx||yyy; xx1={V1 V2}; create pol from d[colNAME=xx1]; append from d ; close pol; xx2={V3 V4}; create perp from kk[colNAME=xx2]; append from kk ; close perp; data pol; set pol; TYPE="pol"; data perp; set perp; TYPE="per"; DATA SSAMMI; SET SUMAS; SSAMMI =COL1; PERCENT =COL2; PERCENAC=COL3; DFAMMI =COL4; DFE =COL5; MSE =COL6; DROP COL1 - COL6; MSAMMI=SSAMMI/DFAMMI; F_AMMI=MSAMMI/MSE; PROBF=1-PROBF(F_AMMI,DFAMMI,DFE); PROC PRINT DATA=SSAMMI NOOBS; VAR SSAMMI PERCENT PERCENAC DFAMMI MSAMMI F_AMMI PROBF; / prepare data for plotting / PROC SORT DATA=RAW; BY GEN; PROC MEANS DATA = RAW NOPRINT; BY GEN ; VAR YLD; OUTPUT OUT = MEDIAG MEAN=YLD; DATA NAMEG; SET MEDIAG; TYPE = 'GEN'; NAME = GEN; KEEP TYPE NAME YLD; PROC SORT DATA=RAW; BY ENV; PROC MEANS DATA = RAW NOPRINT; BY ENV ; VAR YLD; OUTPUT OUT = MEDIAE MEAN=YLD; DATA NAMEE; SET MEDIAE; TYPE = 'ENV'; NAME1 = 'S'||ENV; NAME = COMPRESS(NAME1); KEEP TYPE NAME YLD; DATA NAMETYPE; SET NAMEG NAMEE; DATA BIPLOT0 ; MERGE NAMETYPE SCORES; DIM1=COL1; DIM2=COL2; DIM3=COL3; DROP COL1-COL3; data biplot ; set biplot0 pol perp; PROC PRINT DATA=BIPLOT NOOBS; VAR TYPE NAME YLD DIM1 DIM2 DIM3; Data labels; set biplot ; retain xsys '2' ysys '2' ; length function text $8 ; text = name ; if type = 'GEN' then do; color='black '; size = 1; style = 'hwcgm001'; x = dim1; y = dim2; if dim1 >=0 then position='5'; else position='5'; function = 'LABEL'; output; end; if type = 'ENV' then DO; color='black '; size = 1; style = 'hwcgm001'; x = 0.0; y = 0.0; function='MOVE'; output; x = dim1; y = dim2; function='DRAW' ; output; if dim1 >=0 then position='5'; else position='5'; function='LABEL'; output; end; if type = "per" then do; color='red'; line=2; size =1; style = 'hwcgm001'; x=0.0; y=0.0; function='MOVE'; output; x=v3; y=v4; function='DRAW'; output; end; / graphing the biplot */

Proc gplot data=biplot; Plot dim2dim1 v2v1 / overlay Annotate=labels frame Vref=0.0 Href = 0.0 cvref=black chref=black lvref=3 lhref=3 vaxis=axis2 haxis=axis1 vminor=10 hminor=10 nolegend; symbol1 v=none c=black h=0.7 ; symbol2 v=none c=blue i=j line=5 ; axis2 length = 6 in order = (-1.0 to 1.0 by 0.2) label=(f=hwcgm001 h=1.4 a=90 r=0 'Factor 2') value=(h=0.8) minor=none; axis1 length = 6 in order = (-1.0 to 1.0 by 0.2) label=(f=hwcgm001 h=1.2 'Factor 1') value=(h=0.8) minor=none; run; quit;

DATA RAW;

;

ArminSaed / GGE_SAS

GGE SAS code #1