/*..............uprav.sas.......................*/
/*    priprava souboru uzitkovosti a rodokmenu pro ,odchovny */
/*   upravy s dodatky dalsich efektu                         */
/*........................VUZV  20.12.2003...*/
dm output 'clear';
dm log 'clear';
filename byk 'j:\hanka/b/pol3odchrok/byk.dat';
filename uzitko 'j:\hanka/b/pol3odchrok/uzitkob';
filename rodok 'j:\hanka/b/pol3odchrok/RODOKMEN1';
/*...............cteni uzitkovosti..............*/
data byci; infile byk;
input cislo $ 1-10 odch $ 11-20 otecreg $ 21-30 plemeno $ 31-40 naroz 41-50
dat  51-60  hm 61-70 matka $ 71-80 vek 81-90 _freq_ 91-100 trid 101-110
prirdo 111-120 prirpo 121-130 prupri 131-140 ;
			format dat ddmmyy.;
			format naroz ddmmyy.;
			format rok $4.; format rokvaz $4.;
			format mes $2.; format mesvaz $2.;	
		rok =year(naroz); rokvaz =year(dat);
		mes =month(naroz); mesvaz =month(dat);			proc means;title "uzitkovost";
	data byci; set byci; if otecreg = "" then delete; 
if prirdo > 2300 then delete; if prirpo > 2300 then delete; if prupri > 2300 then delete;
												proc means;title omezeni na uzitkovost;
/*..vahove koeficienty pro vazenou analyzu pro hmotnost a prirustek ..*/
data vah; set byci;
y1 = 5.904137022+0.188527571*vek-0.000420267*vek*vek+0.000000373*vek*vek*vek;
	vahm = 28.37**2/y1**2;
y2 = 248.9253568-0.2473416*vek+0.0006625*vek*vek-0.0000003*vek*vek*vek;
	vap = 228.19**2/y2**2;
												proc means; title "vahove koeficienty";
/*.......koeficienty Leg. polynomu.......................*/
data leg;set vah (drop = y1 y2);
	vmin = 3;
	vmax = 420;
		sv = 2*((vek-vmin)/(vmax-vmin))-1;
p0 = 1;											b0=p0*sqrt(1/2)/sqrt(1/2);
p1 = sv;										b1=p1*sqrt(3/2)/sqrt(1/2);
p2 = 0.5*(3*sv*sv-1);							b2=p2*sqrt(5/2)/sqrt(1/2);
p3 = 0.5*(5*sv**3-3*sv);						b3=p3*sqrt(7/2)/sqrt(1/2);
p4 = 1/8*(35*sv**4-30*sv**2+3);					b4=p4*sqrt(9/2)/sqrt(1/2);
p5 = 1/8*(63*sv**5-70*sv**3+15*sv);				b5=p5*sqrt(11/2)/sqrt(1/2);
p6 = 1/16*(231*sv**6-315*sv**4+105*sv**2-5);	b6=p6*sqrt(13/2)/sqrt(1/2);
												proc means;title"leg polynom";
/*....		precislovani jedincu na plemennou hodnotu..................*/
data prec;set leg (drop = vmax vmin);
	kodje=cislo||otecreg;      
		proc sort;by kodje;
	data cis; set prec (keep =  kodje );jed = kodje;	/*proc means;*/
	data otc; set prec (keep = otecreg );jed = otecreg; /*proc means;*/
data kod; set cis otc; if jed="" then delete; ;	/*proc means;*/
	proc sort; by jed;
			data a;set kod;by jed; if first.jed;				/*proc means;*/
		data a ;set a ;kod = _n_;kodje = jed;keep kod kodje otecreg;
				data aa; set a; drop otecreg;
			proc sort; by kodje;								/*proc means;*/
	/* .............precislovany jedinec do souboru uzitkovosti......*/
data kk;merge prec(in=je) aa;by kodje;if je;	proc means;title"precislovany jedinec"; 
	/*...............precislovany otec do souboru uzitkovosti........*/
proc sort; by otecreg;
	data b; set a; kodot =kod; keep otecreg kodot; if otecreg = "" then delete; 
		proc sort; by otecreg;
	data b; merge kk(in=ot) b; by otecreg; if ot;		proc means;title"precislovany otec";
/* ........ precislovani interakce stanice x plemenik ....................*/
data stan; set b;
	stbyk=odch||kodot;
		proc sort; by stbyk;
        data dat; set stan; by stbyk; if first.stbyk; keep stbyk;
			data dat; set dat; kododby = _n_;
data did;
/*.......................precislovane odch x byk do spouboru uzitkovosti ....*/
	merge stan(in=odb) dat;by stbyk; if odb; 	proc means;title"precislovana odch x otec";
/*........precislovani datumu vazeni.....kodtd.....*/
data vaz; set did; 
	td = dat||odch;
		proc sort; by td;
		data dat; set vaz;by td;if first.td; keep td;
			data dat;set dat; kodtd = _n_;
data dd;
/*.........................precislovany td do souboru uzitkovosti.....*/
	merge vaz(in=va) dat;by td;if va;			proc means; title"precislovany den vazeni"; 
/*......................................vytvoreni efektu rok x odchovna....*/
data stan; set dd;
	odrok=odch||rok;
		proc sort; by odrok; /*proc freq; tables odrok; run; ........................*/
        data dat; set stan; by odrok; if first.odrok; keep odrok;
			data dat; set dat; kododro = _n_;
data did;
/*.......................precislovane odch x rok do souboru uzitkovosti ....*/
	merge stan(in=odb) dat;by odrok; if odb;
										proc means; title"precislovane odcho x rok naroz";
/*.............precislovani roku....cir...........................*/
data ro; set did;
	proc sort; by rok;				
	data ddd; set ro;by rok;if first.rok;keep rok;
		data ddd;set ddd; cir = _n_; 
		/*..................precislovany rok do uzitkovosti..........*/
data uzit;
	merge ro(in=ror) ddd;by rok; if ror;		proc means; title"precislovany rok";
/*.............precislovani odchovny....cod...........................*/
data ro; set uzit;
	proc sort; by odch;				
	data ddd; set ro;by odch;if first.odch;keep odch;
		data ddd;set ddd; cod = _n_; 
		/*..................precislovana odchovna do uzitkovosti..........*/
data uzit;
	merge ro(in=ror) ddd;by odch; if ror;		proc means; title"precislovana odchovna";
/* .......................zapis uzitkovosti .......................*/
data uzitko;set uzit (keep = hm vek prirdo prirpo prupri vahm vap b0 b1
	b2 b3 b4 b5 b6 cislo kod otecreg kodot td kodtd rok cir kododby cod kododro);
			if prirdo < 0 then prirdo = 0.001;
			if prirpo < 0 then prirpo = 0.001;
			if prupri < 0 then prupri = 0.001;
			if prirdo = . then prirdo = 0;
			if prirpo = . then prirpo = 0;
			if prupri = . then prupri = 0;
file uzitko;
put hm    vek  prirdo  prirpo  prupri  vahm 
	vap    b0   b1       b2       b3     b4 
    b5    b6    kod    kodtd   cir    kododby 
	cod   kododro ;
													proc means;title"cela uzitkovost";
												/*	proc freq; tables kododby;*/
/*..........................zapis rodokmenu....................*/
data r1; set uzitko;
	if kodot=. then kodot=0; kot= kodot; kodm=0; keep kod kot kodm;
data r2; set uzitko;
	if kodot=. then delete; kod=kodot; kot=0; kodm=0; keep kod kot kodm;
													proc means;title"otci";			
data rod; set r1 r2; proc sort nodup;by kod;		
data rod; set rod;
file rodok; put kod kot kodm;						proc means;title"rodokmen";
/*............................tisk nekolik radku vysledku......*/
data zko; set uzitko;
if _n_ > 10 then delete;
proc sort; by kod vek;								proc print; 
/*proc export data = uzitk*/
/*outfile = 'uzitkov'*/
/*dbms = tab;*/
/*run;*/
data a; set uzit (keep = hm vek prirdo prirpo prupri vahm vap b0 b1
	b2 b3 b4 b5 b6 cislo kod otecreg kodot td kodtd rok cir kododby cod kododro);
	if prirdo < 0 then prirdo = 0.001;
			if prirpo < 0 then prirpo = 0.001;
			if prupri < 0 then prupri = 0.001;
proc means;title"data pro GLM";
	proc sort; by kodtd;
proc glm;
model hm prirdo prirpo prupri = b0 b1 b2 b3;
proc glm; absorb kodtd;
model hm prirdo prirpo prupri = b0 b1 b2 b3;
proc glm; class cir; absorb kodtd;
model hm prirdo prirpo prupri = cir b1*cir b2*cir b3*cir;
proc glm; class cod;absorb kodtd;
model hm prirdo prirpo prupri = cod b1*cod b2*cod b3*cod;
proc glm; class kododro;absorb kodtd;
model hm prirdo prirpo prupri = kododro b1*kododro b2*kododro b3*kododro;
proc glm; class kododro; absorb kodtd;
model hm prirdo prirpo prupri = kododro b1*kododro b2*kododro b3*kododro kododby;
/*..................................................................*/
proc glm;weight vap;
model prirdo prirpo prupri = b0 b1 b2 b3;
proc glm; absorb kodtd;weight vap;
model prirdo prirpo prupri  = b0 b1 b2 b3;
proc glm; class cir;absorb kodtd;weight vap;
model prirdo prirpo prupri = cir b1*cir b2*cir b3*cir;
proc glm; class cod; absorb kodtd;weight vap;
model prirdo prirpo prupri = cod b1*cod b2*cod b3*cod;
proc glm; class kododro; absorb kodtd;weight vap;
model prirdo prirpo prupri  = kododro b1*kododro b2*kododro b3*kododro;
proc glm; class kododro; absorb kodtd;weight vap;
model prirdo prirpo prupri  = kododro b1*kododro b2*kododro b3*kododro kododby;
/*..................................................................*/
proc glm;weight vahm;
model hm  = b0 b1 b2 b3;
proc glm; absorb kodtd;weight vahm;
model hm  = b0 b1 b2 b3;
proc glm; class cir; absorb kodtd;weight vahm;
model hm = cir b1*cir b2*cir b3*cir;
proc glm; class cod; absorb kodtd;weight vahm;
model hm  = cod b1*cod b2*cod b3*cod;
proc glm; class kododro; absorb kodtd;weight vahm;
model hm  = kododro b1*kododro b2*kododro b3*kododro;
proc glm; class kododro; absorb kodtd;weight vahm;
model hm  = kododro b1*kododro b2*kododro b3*kododro kododby;
run;
/*..............................konec..............................*/