###################################one way random effects model################################################
##install package lme4
library(nlme)
ex.data <- read.table(file="C:/jenn/teaching/stat579/data/influent.txt",header=T)
attach(ex.data)
summary(ex.data)
## influent nitrogen type
## Min. :1.000 Min. : 2.00 Min. :1.000
## 1st Qu.:2.000 1st Qu.:14.00 1st Qu.:1.000
## Median :3.000 Median :20.00 Median :2.000
## Mean :3.162 Mean :20.97 Mean :1.865
## 3rd Qu.:5.000 3rd Qu.:28.00 3rd Qu.:2.000
## Max. :6.000 Max. :42.00 Max. :3.000
nrow(ex.data)
## [1] 37
head(ex.data)
## influent nitrogen type
## 1 1 21 2
## 2 1 27 2
## 3 1 29 2
## 4 1 17 2
## 5 1 19 2
## 6 1 12 2
nrow(ex.data)
## [1] 37
head(ex.data)
## influent nitrogen type
## 1 1 21 2
## 2 1 27 2
## 3 1 29 2
## 4 1 17 2
## 5 1 19 2
## 6 1 12 2
##descriptive statistics
tapply(nitrogen,influent,mean)
## 1 2 3 4 5 6
## 21.55556 13.85714 16.80000 24.50000 14.40000 36.40000
tapply(nitrogen,influent,sd)
## 1 2 3 4 5 6
## 5.746980 7.358183 4.086563 6.348228 9.502631 5.029911
tapply(nitrogen,influent,length) #unbalanced design
## 1 2 3 4 5 6
## 9 7 5 6 5 5
#declare grouping variable as a factor
influent.f<-as.factor(influent)
boxplot(nitrogen~influent.f,main="Boxplot",xlab = "Influent",
ylab = "Nitrogen")
![](data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAABUAAAAPACAMAAADDuCPrAAAAZlBMVEUAAAAAADoAAGYAOmYAOpAAZpAAZrY6AAA6OgA6kNtmAABmOgBmZjpmZmZmtrZmtv+QOgCQkDqQkNuQ29uQ2/+2ZgC2kDq225C2/7a2///bkDrb/9vb////tmb/25D//7b//9v///+RQCuMAAAACXBIWXMAAB2HAAAdhwGP5fFlAAAgAElEQVR4nO3d62LjNqOeUSdNJ530MGkbdzdunMP932Ql0QdapicUTBAvwLX+7P1lxiIIkM/oQNF3/wBQ5K71AAB6JaAAhQQUoJCAAhQSUIBCAgpQSEABCgkoQCEBBSgkoACFBBSgkIACFBJQgEICClBIQAEKCShAIQEFKCSgAIUEFKCQgAIUElCAQgIKUEhAAQoJKEAhAQUoJKAAhQQUoJCAAhQSUIBCAgpQSEABCgkoQCEBBSgkoACFBBSgkIACFBJQgEICClBIQAEKCShAIQEFKCSgAIUEFKCQgAIUElCAQgIKUEhAAQoJKEAhAQUoJKAAhQQUoJCAAhQSUIBCAgpQSEABCgkoQCEBBSgkoACFBBSgkIACFBJQgEICClBIQAEKCShAIQEFKCSgtPXnz3dvfdvy0f/65fSIP/z2b3/t8euWG+U4BJS23gV004SuCuhpCAJKEQGlrYWAbpizFQH9+9dNt8ihCChtLQX0319zr7UioPfbJptDEVDaugT0NXGXnN192erRBZSqBJS2rgI69ezH3zd6dAGlKgGlreuAXv/vzxFQqhJQ2loM6OwZ6OUzntl/uSTx+QceXj5zenj6+P7h7ef41wG9frTHWtdPcQwCSlvXAb0k7ac/5n/69pOlqXmXvzH96eX/nQL60sPnQr4N6PtHE1A+RUBpa/E90OdX1A9LVzfdP/+P6enk9LOXv/nf332Q/yagC48moHyKgNLW24BOTXx+/vh4d/e+ea9/53H2Hx+u/urTY8wDuvRoAsqnCChtLVwH+ub19zx0TyWcfuTL9MdPr/YfXp94zn5sHtCPHs2HSJQTUNp6H9CXi0DvZ5mb8vf0R7Onm29fnk81nT2NnQX0o0cTUMoJKG0tPAN9U7nntD3O/uTps/T5H7/59H16zG//zAP64aMJKOUElLaWvso5+4Dn5YKjKZpfr37o5eP6h/lffv27rwH98NEElHICSltXn8I/fapzfvb48KaQU+ieX97P3vKc/YeXv3z/PqAfPpqAUk5AaWvxQvpL6R7mxbz6n08v4r8u/+lrHq8CuvRoAko5AaWtd1/dfHh+avm9gD69hr96Biqg7ExAaetdQB+fX8N/5yX8y6dIb98D9RKenQkobX0c0KWPfZ4+Zn+9junNp/C3foh0eTQBpZyA0tbyS/hz2z6+jGnhG/LvL2N6vXjJZUzUIqC0dR3QKXSXp4ofXfr+dKH8f7zeSuQ5oNMzzOnPL3/gQnqqElDaugro05PLL7P//92XL5+fbb7ezO7NdU3zT+hnAf3o0e7fvjkKNxBQ2lq8kP5NKV99m/3El3/e341pbkri/GYiy4/2+p89C+VmAkpbiwF9vjHS0u3sZt90v74f6H/6r9cJ/rfb2c0HIKDcTEBpaymg35b+9DmD97PaXd2R/qf/93x509obKr8+4Ia/yY7jEFDaehfQq99gtPxLOJ7es3z6w9lFo1dvqf7br/SYDWG736XMcQgoQ7i6Th52IaAMQUBpQUAZgoDSgoAyBAGlBQFlCAJKCwLKEASUFgSUIQgoLQgoQxBQWhBQgEICClBIQAEKCShAIQEFKCSgAIUEFKCQgAIUElCAQgIKUEhAAQoJKEAhAQUoJKAAhQQUoJCAAhQSUIBCAgpQSEABCgkoQCEBBSgkoACFBBSgkIACFBJQgEICClBIQAEKCShAIQEFKCSgAIUEFKCQgAIUElCAQgIKUEhAAQoJKEAhAQUoJKAAhQQUoJCAAhQSUIBCAgpQSEABCgkoQCEBBSgkoACFBBSgkIACFBJQgEICClBIQAEKCShAIQEFKCSgAIUEFKBQdkDvADazfaI2f8QNtZ5tYCybN2rrB9xShX8wgMMSUIBCAgpQSEABCgkoQCEBBSgkoACFBBSgkIACFBJQgEICClBIQAEKCShAIQEFKCSgAIUEFKCQgAIUElCAQgIKUEhAAQoJKJAr8he/zYe3+SNu/YBbElDoySb9FNCtCCiMpe05LaBAxwR0RwIKYxHQHQkojEVAdySgMBYB3ZGAwlgEdEcCCmMR0B0JKIxFQHckoMB2BBSgUH8BvX/5ftaX239YQIHtdBbQh6vvuN7aUAEFttNVQP/8+d1dAn78/aZHEFDIsc2tQjbwiT3YcDqmR9z6AV/8/evdm2L+9cvpf//0xy0PIaAQo3E158p3YcsJuTzi1g/44vFdLs9J/XbLQwgoxIg5HY8R0Pv3L9hPT0Jvehs0ZsWAmNPxEAE9Pd38+u4/Ptz2Gj5mxYCY0/EQAT0923z/cv3xto+RYlYMiDkdBXSlmBUDYk7HQwTUS3gYSszpeIiA+hAJhhJzOh4joMuXMb1/VvodMSsGxJyOxwjo4oX0t30VKWbFgJjT8RgBnYr51g+/3fQIMSsGxJyOBwno/FZMEzcTgW7FnI6HCeg/bmcHo4g5HY8U0BtsedMAYGMxp6OALhJQCBZzOh4xoH/9ctsVTBcxKwbEnI4CulLMigExp6OArhSzYkDM6XiIgC5cBXp365X0MSsGxJyOArp2cCkrBsScjocI6Pm78AIKw4g5HY8R0Mtz0JfbiXgPFPoWczoeJKCXryE9f/tdQKFvMafjYQJ6+dXwUzcFFPoWczoeJ6CXe9pdXsYLKPQt5nQ8UEDPv8Xj8svgBRT6FnM6Hiqgl5fxXwQUOhdzOh4roJeX8T/+XwGFrsWcjgcL6PQy/sZfh3QRs2JAzOl4uIBeXsYLKPQs5nQ8XkAvl4QKKHQs5nQ8YkCLxKwYEHM6CuhKMSsGxJyOArpSzIoBMaejgK4Us2JAzOkooCvFrBgQczoK6EoxKwbEnI4CulLMigExp6OArhSzYkDM6SigK8WsGBBzOgroSjErBsScjgK6UsyKATGno4CuFLNiQMzpKKArxawYEHM6CuhKMSsGxJyOArpSzIoBMaejgK4Us2JAzOkooCvFrBgQczoK6EoxKwbEnI4CulLMigExp6OArhSzYkDM6SigK8WsGBBzOgroSjErBsScjgK6UsyKATGno4CuFLNiQMzpKKArxawYEHM6CuhKMSsGxJyOArpSzIoBMaejgK4Us2JAzOkooCvFrBgQczoK6EoxKwbEnI4CulLMigExp6OArhSzYkDM6SigK8WsGBBzOgroSjErBsScjgK6UsyKATGno4CuFLNiQMzpKKArxawYEHM6CuhKMSsGxJyOArpSzIoBMaejgK4Us2JAzOkooCvFrBgQczoK6EoxKwbEnI4CulLMigExp6OArhSzYkDM6SigK8WsGBBzOgroSjErBsScjgK6UsyKATGno4CuFLNiQMzpKKArxawYEHM6CuhKMSsGxJyOArpSzIoBMaejgK4Us2LAXY7yXdhyQi6PuPUDbklAIUbras6U78KWE3J5xK0fcEsCCjFaV3OmfBe2nJDLI279gFsSUIjRupoz5buw5YRcHnHrB9ySgEKM1tWcKd+FLSfk8ohbP+CWBBRitK7mTPkubDkhl0fc+gG3JKAQI+Z0FNCVYlYMiDkdBXSlmBUDYk5HAV0pZsWAmNNRQFeKWTEg5nQU0JViVgyIOR0FdKWYFQNiTkcBXSlmxYCY01FAV4pZMSDmdBTQlWJWDIg5HQV0pZgVA1p/f3OmfBe2nJDLI279gFsSUIjRupoz5buw5YRcHnHrB9ySgMJY2p7TAgp0TEB3JKAwFgHdkYDCWAR0RwIKYxHQHQkojEVAdySgMBYB3ZGAwlgEdEcCCmxHQAEKCShAIQEFKCSgAIUEFKCQgAIUElCgY64D3ZGAwlgEdEcCCmMR0B0JKIxFQHckoDAWAd2RgMJYBHRHAgpjEdAdCSiMRUB3JKAwFgHdkYAC2xFQgEICClBIQAEKCShAIQEFKCSgAIUEFOiY60B3JKAwFgHdkYDCWAR0RwIKYxHQHQkojEVAdySgMBYB3ZGAwlgE9Db3d8++3P7DAgpjEdAbPNy9dWtDBRTGIqCr/fnz3bUff7/pEQSUet4dnQVa7wO36Smgf//6tph//XL63z/9cctDOECpZot+Oj4701NAH9/l8pzUb7c8hAOUdhx94+kpoPfvX7CfnoTe9DaoQ5h2HH3j6Sigp6ebX9/9x4fbXsM7hGnH0TeejgJ6erb5/uX6420fIzmEacfRNx4BhZ04+sbTUUC9hKdvjr4aXAe6lg+R6JqjrwYBXWv5Mqb3z0q/wyFMO46+GgR0rcUL6W/7KpJDGMYioKtdivnWD7/d9AgCCmMR0BvcX/XTzUTg2AT0Nm5nB7wQ0GrcrIHNbHKrkA20noe9hU+agMIKm5zHm2g9E/tKn7NeA3r5RP62S5jOjnb4sZWUIydlHEw6C+jj3XQt6OPTvyw33k/Z4UehlCMnZRxMugro9LTz7offXm9N74707CLlyEkZB5OuAvr8AfxP/+epnPfuSM8+Uo6clHEw6Smgp+ed52xefrHcUzcf3ZGeXaQcOSnjYNJTQB+eXrCfnnc+fwHp9KLezUTYQcqRkzIOJh0F9OV2dqdnoi8v3N3Ojl2kHDkp42DSUUBfbqg8f9rphsrsIuXISRkHkx4DenoNL6DsK+XISRkHk74C+nTh/P38JbyAsoOUIydlHEw6Cujplfv79zvvfYjEHlKOnJRxMOkooOfrl66vWXp0R3p2kXLkpIyDSU8BPd9P+cv1f0i6I332bQ/4jJSVSRkHk54COn0D/vVJ6P3b/7lG1cNvi346PUKlrEzKOJh0FdBzQWe/w+Ph5n4Od/NV9pKydinjYNJXQE/POmcv2R9uvye9gFImZe1SxsGkt4B+koBSJmXtUsbBREAPs3U+I2XtUsbBREAPs3U+I2XtUsbBREAPs3U+I2XtUsbBREAPs3U+I2XtUsbBREAPs3U+I2XtUsbBREBhhZQjJ2UcTAQUVkg5clLGwURAYYWUIydlHEwEFFZIOXJSxsFEQGGFlCNno3G4781GBBRWSDlythnHFv0MmZDGBBRWSDlydhpHyu7GE9DDbJ3PSFk7Ac0ioIfZOp+RsnYCmkVAD7N1PiNl7QQ0i4AeZut8RsraCWgWAT3M1vmMlLUT0CwCepit8xkpayegWQT0MFvnM1LWTkCzCOhhts5npKxdyjiYCOhhts5npKxdyjiYCCiskHLkpIyDiYDCCilHTso4mAjo+h/NsOFssF7KxKeMg4mArv7JFFvOB2ulzHvKOJgIaPWf3FbKOI4mZd5TxsFEQKv/5LZSxnE0KfPuOtAsAlr9J7eVMo6jSZl3Ac0ioNV/clsp4zialHkX0CwCWv0nt5UyjqNJmXcBzSKg1X9yWynjOJqUeRfQLAJa/Se3lTKOo0mZdwHNIqDVf3JbKeM4mpR5F9AsAlr9J7eVMo6jSZl3Ac0ioNV/clsp4zialHkX0Cw7BfTPn+dfRvzx9603upaAUiZl3lPGwWSfgD68/Ta3gJZLGcfRpMx7yjiY7BLQxzsB3UrKOI4mZd5TxsFkl4Det2zmGwJKmZR5TxkHkz0C+vevKf0UUAqlzHvKOJjsEdC/frn7uvVWCgkoZVLmPWUcTAS0+k9uK2UcR5My7ynjYLLTS3gB3UrKOI6m6i8ZuMlOu7vLZvq3y4dID3c//bH1ZsoIKGVaZ/PVTru7y2b6t0tAc17DCyhlWmfz1U67u8tm+rfPhfR//hzyHFRAKdM6m6922t1dNtO/nT5EesuF9OVSxnE0TVq5aKfd3WUz/RPQ6j+5rZRxHE2TVi7aaXd32Uz/BLT6T24rZRxHkzLvAprF7eyq/+S2UsZxNCnzLqBZBLT6T24rZRxHkzLvAppFQKv/5LZSxnE0KfOeMg4muwX0/EboD7/989D2ciYBpUzKvKeMg8lOAZ0+RzoH9O7u29ZbvIGAUiZl3lPGwWSfgD59Dn8K6P1d04IKKGVS5j1lHEx2Cejfv97d/fTHnz+fAtr45soCSpmUeU8ZB5O9fqXH1/P3Oc8BPb+Ib/fFeAGlTMq8p4yDyU6/0uP80dFTQE9PR9t9kCSglEmZ95RxMNnxfqBPAT09BfVNpHIp4zialHl3HWiWnb7Kef7c6DmgjwL6CSnjOJqUeRfQLAJa/Se3lTKOo0mZdwHN0uAl/L33QD8hZRxHkzLvApplpw+Rzs85nwJ6ej76ZettriaglEmZdwHNstdlTF+eA3q+JrTdlfQCSpmUeRfQLPt8E+ly9fwloA93dy1/u4eAUiZl3gU0y65f5Wx9O2UBpVTKvAtolp1uJnJ+4R7QTwGlUMq8C2iW3W5n95TQpvkUUEqlzHvKOJi4oXL1n9xWyjiOJmXeU8bBRECr/+S2UsZxNCnznjIOJgJa/Se3lTKOo0mZ95RxMBHQ6j+5rZRxHE3KvKeMg8lO34W/1upSUAGlTMq8p4yDSaOA3jX6OqeAUiZl3lPGwaRZQNsUVEApkzLvrgPNss97oH/+fDfdiekS06/z/7AvAaVMyrwLaJZdAnrK5euvQXq8pLPRPZkElDIp8y6gWXa6nd2Xd//rockHSQJKmZR5F9AsO96R/tl0R/o296UXUMqkzLuAZhHQ6j+5rY3GsfSx3q02GUgvUnZXQLPs9Cs93r6EP794F9Ay24xji36GTMhOUnZXQLPs9B7omw+RLjlt85uRBDRrM/1ImRDrn2WvT+Gfn2+eL2M6fwr/MG/qfgQ0azP9SJkQ659ln+tAH9689Pt2KarrQIs4gdpImZCUcTDZ6WYi52I+OYfz9DK+yRNQAQ3bTD9SJiRlHEx2uxvT0/c5p1fyj22efwpo2mb6kTIhKeNg4nZ21X9yWwLaRsqEpIyDiYBW/8ltCWgbKROSMg4mbV7CNyOgWZvpR8qEpIyDyU4Bnd3Rzq81/pSUcRxNyrynjIPJPgF9c0fQRp8fXQgoZVLm3SuQLLsE9PI74aevw58vYGr1+zz+EVBKpcy7gGbZJaCPs2ieY/rt/V/ZiYBSJmXeBTTLTt+Fn73x2ehWyhMBpUzKvAtolp3uxjT/2lGbWylPBJQyKfMuoFma3Q+0DQGlTMq8C2gWAa3+k9tyArWRMiHWP4uX8NV/cltOoDZSJsT6Z/EhUvWf3JYTqI2UCbH+WVzGVP0nt+UEaiNlQlLGwcSF9NV/clsC2kbKhKSMg0mDr3K2/DK8gGZtph8pE5IyDiZuJlL9J7cloG2kTEjKOJi4nV31n9yWgLaRMiEp42DihsrVf3JbAtpGyoSkjIPJPpcxbfm88/7lzYCCq6EElDIp854yDiY7XUi/1Qfvb38/8u0NFVDKpMy7VyBZdvoq5zaXzs9+OXLhR1ICSpl3R14zO+3uLpvpX0/PQC+Xk779TtOtF5UKKGXaRnNup93dZTP92+ubSFt8+ejxXS5v/lqTgNJOR2vX0VDb2udT+NNr7w2eg96/f8F+67sDAko7Ha1dR0Nta6f3QN8q+1D+6qZOkxtv7SSgtNPR2nU01LY6CujVbUUnN95cVECzNkMo67+SgFb/yW0JKDuw/it19E0kL+HPBJQdWP+VOgqoD5HOBBRy9BTQ5cuY3j8r/Q4BzdoM9G2f90B/+G32P+9LL2lavJD+tvdTBTRrM9C3JgEtvbfIu0+j7u7ePPKKwQlo1Gagb/sH9FTB8psz3V/1c8+biaQo3YMbd3eXzRyLSR1P5YBe3z1p8qlvJTW7nV2K0j24cXd32cyxmNTxVA7owovuuxs/9/nMUDbMT90o3mLLGWqt9Vy+2Glvd9kMO6r9Ev5x4WDdq58Cmq71VM7ss7t7bGUTHQ21rQYfIrUkoFFidkZAr3Q01LY6C+jllsqzt1D//nW/r3KmKN2DQDE7I6BXOhpqWz1dSP/6AdLLx0d7BrT0J7eVMo5NxOyMgF7paKhtdRXQ1w/gn5+ECmjXYnZGQK90NNS26gb0/CH8t63uxnR5/X6+H9PDa0EFtGsxOyOgVzoaals9BfT++XtH5+90TgUV0OW/k2GTndlFzEBSmJCVOgroKZYv733ePxVUQBf/SooNdmYfMQNJYUJW6ug90Dc3VH4qqICW/ZVdCGi/TMhKvQb0XNAvAlr6V3YhoAyv24BOtwIV0LK/sgsBZXgdBXT+Huj0P+++CWjZX9nFmoDG2GNCGE9HAT1fvjT/rXJ//nz3w/8W0KK/sgsBZXj1P4Vf8InrQOc3Inm8/cEO8VXOlBoI6M0TQm96Cuh0d9HZy/jHmx9MQHckoDdPCL3pKqCXgs6fg56fk+4U0C0c7LaTAnrzhNCbnd8DvdxN6RM3BH33u+EfBLTVZv6VgN48ITE6Gmpbuwb08ms1P/cLPT5JQHfkMqYmW9lER0Nta8+AXt7CbHtv5UMcFyk7KaBNtrKJjoba1n4Bnd4PLfhFcFs6xHGRspMC2mQrm+hoqG3tFtDLvTzLf5/xRg5xXKTspIA22comOhpqWzsFdLpk89vSH+3qEMdFyk4KaJOtbKKjoba1S0AjXr1fHOK4SNlJAe2XCVlpj4BmvHq/OMRxkbKTAtovE7JS/YB+9tLPTR3iuEjZSQHtlwlZqXZA21/6+YbrQHckoAyvckADLv18Q0B3JKAMr6/vwn9+cAK6HwFtpNbXXW/Veh72IKDjbT3lyBXQNnbO5He0nokdCOh4W085cAW0yVZiZjVlHFXtfDem1gR0RwLaZCsxs5oyjqoEdLytpxy4AtpkKzGzmjKOqgR0vK2nHLgC2mQrMbOaMo6qBHS8raccuALaZCsxs5oyjqoEdDwpOymgTbYSM6sp46hKQMeTspMC2mQrMbOaMo6qBHQ8KTspoE22EjOrKeOoSkDHk7KTAtpkKzGzmjKOqgR0PCk7KaBtpOxMyjiqEtDxpOykgLaRsjMp46hKQMeTspMC2kbKzqSMoyoBHW/rKQeugLaRsjMp46hKQMfbesqBK6BtpOxMyjiqEtDxtp5y4ApoGyk7kzKOqgR0vK2nHLgC2kbKzqSMoyoBHW/rKQeugDbZSsyspoyjKgEdb+spB66ANtlKzKymjKMqAR1v6ykHroA22UrMrKaMoyoBHW/rKQeugDbZSsyspoyjKgEdb+spB66ANtlKzKymjKMqAR1Pyk4KaJOtxMxqyjiqEtDxpOykgDbZSsyspoyjKgEdT8pOCmiTrcTMaso4qhLQ8aTspIA22UrMrKaMoyoBHU/KTgpoGyk7kzKOqgR0PCk7KaBtpOxMyjiqEtDxpOykgLaRsjMp46hKQMfbesqBK6BtpOxMyjiqEtDxtp5y4ApoGyk7kzKOqgR0vK2nHLgC2kbKzqSMoyoBHW/rKQeugLaRsjMp46hKQMfbesqBK6BNthIzqynjqEpAx9t6yoEroE22EjOrKeOoSkDH23rKgSugTbYSM6sp46hKQMfbesqBK6BNthIzqynjqEpAx9t6yoEroE22EjOrKeOoSkDHk7KTAtpkKzGzmjKOqgR0PCk7KaBNthIzqynjqEpAx5OykwLaZCunzaTYZXfbEtDxpOykgDbZioDuSUDHk7KTAtpG62y+aj0TOxDQ8aTspIC20Tqbr1rPxA4EdDwpOymgbbTO5qvWM7EDAR1v6ykHroC20Tqbr1rPxA4EdLytpxy4AtpGys6kjKMqAR1v6ykHroC2kbIzKeOoSkDH23rKgSugbaTsTMo4qhLQ8baecuAKaJOtxMxqyjiqEtDxtp5y4Apok63EzGrKOKoS0PG2nnLgCmiTrcTMaso4qhLQ8baecuAKaJOtxMxqyjiqEtDxtp5y4Apok63EzGrKOKoS0PGk7KSANtlKzKymjKMqAR1Pyk4KaJOtxMxqyjiqEtDxpOzkmoDG2GaHI8bRz/oPQEDHk7KTRwtozjg2eJANpIyjKgEdT8pOpozjaFLmPWUcVQnoeFJ28mBv+cVImZCUcVQloONJ2UkBbSNlQlLGUZWAjrf1lAP3YLsbI2VCUsZRlYCOt/WUA/dguxsjZUJSxlGVgI639ZQD92C7GyNlQlLGUZWAjrf1lAP3YLsbI2VCUsZRlYCOt/WUA/dguxsjZUJSxlGVgI639ZQD92C7GyNlQlLGUZWAjrf1lAP3YLsbI2VCUsZRlYCOt/VNvlK4iV12lysp854yjqoEdLytt87mq112lysp854yjqoEdDyts/mq9UwcU8q8p4yjKgEdT+tsvmo9E8eUMu8p46hKQMfTOpuvWs/EMaXMe8o4qhLQ8bTO5qvWM3FMKfOeMo6qBHQ8KTuZMo6jSZn3lHFUJaDjSdlJ14G2kTIhKeOoSkDHk7KTAtpGyoSkjKMqAR1v6ykH7sF2N0bKhKSMoyoBHW/rKQfuwXY3RsqEpIyjKgEdb+spB+7BdjdGyoSkjKMqAR1v6ykH7sF2N0bKhKSMoyoBHW/rKQfuwXY3RsqEpIyjKgEdb+spB+7BdjdGyoSkjKMqAR1v6ykH7sF2N0bKhKSMoyoB3fTRI779mHLgpozjaFLmPWUcVQnolg8uoHMp4zialHlPGUdVAjqelJ1MGcfRpMx7yjiqEtDxpOxkyjiOJmXeU8ZRlYCOJ2UnU8ZxNCnznjKOqgR0PCk7mTKOo0mZ95RxVCWg40nZyZRxHE3KvKeMoyoBHU/KTroOtI2UCUkZR1UCOp6UnRTQNlImJGUcVQnoeFJ2UkDbSJmQlHFUJaDjSdlJAW0jZUJSxlGVgI4nZScFtI2UCUkZR1X9BfT+5VuPX27/YWu6IwFtI2VCUsZRVWcBfbj65vitDbWmOxLQNlImJGUcVXUV0D9/fnfvjR9/v+kRrOmOthlHyi1a+pGyuynjqKqngP7969ti/vXL6X//9MctD2FNd7TRraUE9EYpu5syjqp6Cujju1yek/rtloewpjtKGcfRpMx7yjiq6img9+9fsJ+ehN70Nqg13VHKOI4mZd5TxlFVRwE9Pd38+u4/Ptz2Gt6a7ihlHEeTMu8p46iqo4Cenm2+f7n+eNvHSNZ0RynjOJqUeU8ZR1UCOp6UnUwZx9GkzHvKOKrqKKBewq+UspMp4zialLbig7sAAA1OSURBVHlPGUdVHQXUh0grpexkyjiOJmXeU8ZRVU8BXb6M6f2z0u+wpjtKGcfRpMx7yjiq6imgixfS3/ZVJGu6o5RxHE3KvKeMo6qeAjoV860ffrvpEazpjlLGcTQp854yjqq6Cuj8VkwTNxNZkLKTKeM4mpR5TxlHVZ0F9B+3s/t3KTuZMo6jSZn3lHFU1V9Ab3DQ20psc/+NLbSeiWNKmfeUcVQloOPZO5Mfaz0Tx5Qy7ynjqGrogL53iDXl2FIO8pRxVCWgMJaUgzxlHFUJKIwl5SBPGUdVAgpjSTnIU8ZRlYDCWFIO8pRxVNVRQBe+h3TmdnYwl3KQp4yjKgGFsaQc5CnjqKqjgC79VmMBrcRE9Stl7VLGUVVPAb08B73p7nXvHGJNt2Ci+pWydinjqKqrgJ4LeuPtl64cYk23YKL6lbJ2KeOoqq+Anl/F3/QrPK4dYk23YKL6lbJ2KeOoqrOA/vNwt/Cb5dY7xJpuwUT1K2XtUsZRVW8BPb2I/8xT0EOs6RZMVL9S1i5lHFX1FtBPPgU9xJpuwUT1K2XtUsZRVXcB/ZxDrOkWTFS/UtYuZRxVCShLTFS/UtYuZRxVCSiMJeUgTxlHVQIKY0k5yFPGUZWAwlhSDvKUcVQloDCWlIM8ZRxVCSiMJeUgTxlHVQIKY0k5yFPGUZWAwlhSDvKUcVQloCwxUf1KWbuUcVQloCwxUf1KWbuUcVQloCwxUf1KWbuUcVQloCwxUf1KWbuUcVQloCwxUf1KWbuUcVQloCwxUf1KWbuUcVQloCwxUf1KWbuUcVQloCwxUf1KWbuUcVQloCwxUf1KWbuUcVQloDCWlIM8ZRxVCSiMJeUgTxlHVQIKY0k5yFPGUZWAwlhSDvKUcVQloDCWlIM8ZRxVCSiMJeUgTxlHVQIKY0k5yFPGUZWAssRE9Stl7VLGUZWAssRE9esuRuuZ2IGAssRE9at1Nl+1nokdCChLTBSsIKAsMVHHZv1XElCWmKhjs/4rCShLTNSxWf+VBJQlJurYrP9KAsoSE3Vs1n8lAQWuOVFWElDgmhNlJQEFrjlRVhJQ4JoTZSUBBSgkoACFBBSgkICyxETBCgLKEhMFKwgoS0wUrCCgLDFRsIKAssREHZv1X0lAWWKijs36rySgLDFRx2b9VxJQlpioY7P+KwkoS0zUsVn/lQQUuOZEWUlAgWtOlJUEFLjmRFlJQIFrTpSVBBSgkIACFBJQgEICyhITBSsIKEtMFKwgoCwxUbCCgLLERMEKAsoSE3Vs1n8lAWWJiTo267+SgLLERB2b9V9JQFlioo7N+q8koCwxUcdm/VcSUOCaE2UlAQWuOVFWElA4nrsNtN6HCAIKh7NFP51KZwIKUEhAAQoJKEAhAWWJiYIVBJQlJgpWEFCWmChYQUBZYqJgBQFliYmCFQSUJSYKVhBQlpgoWEFAWWKiYAUBZYmJghUEFKCQgAIUElCAQgIKUEhAAQoJKEAhAQUoJKAsMVGwgoCyxETBCgLKEhMFKwgoS0wUrCCgLDFRsIKAssREwQoCyhITBSsI6DHdbaD1PkBzAnpIW/TTVIKAAhQSUIBCAgpQSEABCgkoQCEBBSgkoACFBBSgkIACFBJQgEICClBIQAEKCShAIQEFKCSgAIUEFKCQgAIUElCAQgIKUEhAAQoJKEChwwUUYDubN2rrB9xS68kGxrJ5o7Z+wHF09AaCoVbQz0gNtZ3BdmdLHS21oVbQz0gNtZ3BdmdLHS21oVbQz0gNtZ3BdmdLHS21oVbQz0gNtZ3BdmdLHS21oVbQz0gNtZ3BdmdLHS21oVbQz0gNtZ3BdmdLHS21oVbQz0gNtZ3BdmdLHS21oVbQz0gNtZ3BdmdLHS21oVbQz0gNtZ3BdmdLHS21oVbQz0gNtZ3BdmdLHS21oVbQz0gNtZ3BdmdLHS21oVbQz0gNtZ3BdmdLHS21oVbQz0gNtZ3BdmdLHS21oVbQz0gNtZ3BdmdLHS21oVbQz0gNtZ3BdmdLHS21oVbQz0gNtZ3BdmdLHS21oVbQz0gNtZ3BdmdLHS21oVbQz0gNtZ3BdmdLHS21oVbQz0gNtZ3BdgdgPwIKUEhAAQoJKEAhAQUoJKAAhQQUoJCAAhQSUIBCAgpQSEABCgkoQCEBBSgkoACFBBSgkIACFBJQgEICClBIQAEKCShAIQEFKCSgAIUEFKCQgAIUElCAQgIKUEhAAQoJ6Mf++uXH31uPYYU/f747+eG31uNYYRrq19bDWOs03G+tx/Bv/v717kX8YP957GSc6wnoh06HZgcB/euX59MnfrCv53ofZ9B5auNH+rr++YOd/v3s4VBdT0A/dN/DQs/Pn/DR9vVc6Z/LAZA/0Me7bmb1pZ93dz/90XowWxHQD1xO9/AknTw8vyQ+n+1fWo/mu56Hem5+DyfQYwdNusxq/BgnL+v+csyOQECXTf9axgd09irzIXy4p3+Rnt6nPQ26g3dsp+f28XG6z171mYeXfzcf+/gXdBUBXXI5ef5zB++BPr4+7Tw/ZU4+3U//JD0PtYdnTed3wP9n/jhPw+ykRa//gM7/3+4J6JLzi4xvPXyI9DB7MXSff7o/6WGk58h3EPrTP/XZb9y8mP0DOhIBXfJwXuseAjp338s/66eTPn5iT6f71x6eKZ9egny7vF6KX/vHgd74nBHQD3UW0NNp1MeLuYcO3lqcXhl3ENCHux/+Sx+fbF8mc7oSI31SbyGgH+osoA9d/At/38P1Ns/P5jsI6P3sKqbsgp6n9KGPC+5uIaAf6iugp9ec2SfQpI+APr3czA/o+Rnd9OL9Mf3aoFNA/1svlyzfQEA/1FVAT/3sYbDPF9Nnt/753ZD8gM6uCAs/AqaVv0zrY/wlyzcQ0A/1FNDHnv5RP3/oEX0CPV9bmR/QuexrGy4BfXqKfGp9/GdeawnohzoKaFf9jL+S/qWbfQU0/E3w+9nrjuzW30RAP9RPQO/TXxNfiz7XXy9Y7Cug4dcJ3c+GFz7UWwjoh3oJ6PnFUfQr4veiT6CHu7dyR/pW9KS+/UczfKi3ENAPdRLQ87f2ezsao0+gXgMa/nx59q3j7PW/jYB+qI+AdvOGfDffOu0poLNJTf+C+fwqgej1v42AfqiLgIZ/HjMzO4HCr7h5Ef6c7p83M/mY/kbO/csAe1n/NQT0Qz0ENP0WTDPnob5cHRT9tO5FfkDfTGr4v6Tn95pergPtYv3XENAP9RDQx7evNqNP99kNycOfKz3LD+ibSU2P0uy9kb4uGvkeAf1QBwGd/5aM+IDOfv1IH/3sIaCzgsYP9fWf+3H6KaAf6yCg89+I1MM5NJ3t3Zw+PQT0OUt9/Jt0uRVC+ll1EwEFKCSgAIUEFKCQgAIUElCAQgIKUEhAAQoJKEAhAQUoJKAAhQQUoJCAAhQSUIBCAgpQSEABCgkoQCEBBSgkoACFBBSgkIACFBJQgEICClBIQAEKCShAIQEFKCSgAIUEFKCQgAIUElCAQgIKUEhAAQoJKEAhAQUoJKAAhQQUoJCAAhQSUIBCAgpQSEABCgkoQCEBBSgkoACFBBSgkIACFBJQgEICClBIQAEKCShAIQEFKCSgAIUEFKCQgAIUElCAQgIKUEhASfPXL3c//fHhnz7cnfz0/369+/H3T2zif5X/LLwSUNJ8L6B//3ru592Xvz8V0Pu7r8U/CzMCSprvBfTx0s+7r58K6OlBBJRNCChpvhfQ0wv4y58JKBEElDT/EtBv5/8roEQQUNIIKN0QUNI8B3SK5PSp+x9Pf/D0FuhzQE//54ffpp/68+eX7E5vlH6b/se7R3l6G/X5z+ETBJQ084D+x8+z3K0L6J/PPzI9RX33KALKdgSUNLOAvjp3clVAX/r5VNB3jyKgbEdASfM2oOfOPdw9v2t59R7o+4Cef+jyn86d/LL8KN4DZSsCSpp5QJ8+KXqYYrgioI/Pb3WeH+b8ZwuPIqBsRUBJMw/oU+he3uD814Dev/yHp04uPIqAshUBJc0soO8+Ifq3gM4vgTr9ly+LjyKgbEVASfPJgM6d/ouAUpGAkuYzAZ19Bi+g1CegpPlkQK++xSSgVCSgpPnkS/irr3gKKBUJKGluCOjsQ/eH5xfsV1fICygVCShpbgvo0wWi50+Pzn/l4eWqz9fLmASUWgSUNLcE9PnbRY9PHxlNH8Nf/sqprU8X0i8E9Du/MwTWE1DS3BLQ18uW/scvL3V88XX5UaaP6n0Xns8TUNLcEtCXXn59uYT+9Uqmr/8sP8r09Xiv4vk8ASXNTQG9vFQ//73Zd5AuN/98em908VEuBf2yz94wNAEFKCSgAIUEFKCQgAIUElCAQgIKUEhAAQoJKEAhAQUoJKAAhQQUoJCAAhQSUIBCAgpQSEABCgkoQCEBBSgkoACFBBSgkIACFBJQgEICClBIQAEKCShAIQEFKCSgAIUEFKCQgAIUElCAQgIKUEhAAQoJKEAhAQUoJKAAhQQUoJCAAhQSUIBCAgpQSEABCgkoQCEBBSgkoACFBBSgkIACFBJQgEICClBIQAEKCShAIQEFKPT/Aeyo3MFh6FGGAAAAAElFTkSuQmCC)
# fit one-way ANOVA, assuming only random error
fit1 <- aov(nitrogen~influent.f)
summary(fit1)
## Df Sum Sq Mean Sq F value Pr(>F)
## influent.f 5 1925 385.0 9.044 2.2e-05 ***
## Residuals 31 1320 42.6
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
plot(fit1$resid~fit1$fitted)
![](data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAABUAAAAPACAMAAADDuCPrAAAAUVBMVEUAAAAAADoAAGYAOpAAZmYAZpAAZrY6AAA6OgA6kJA6kNtmAABmtv+QOgCQkNuQ2/+2ZgC2///bkDrb/7bb/9vb////tmb/25D//7b//9v///+8iDO7AAAACXBIWXMAAB2HAAAdhwGP5fFlAAAgAElEQVR4nO3dbYOaSniAYZO0SZtsGmMb8/L/f2hFRXFlWXgcYWa4ri9nsweMhOVeXgbc/AMgZLP0GwAolYACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKEBQ3gHdACSTPlHJXzGhpf+1gbokb1TqF0zpCb8wgNUSUIAgAQUIElCAIAEFCBJQgCABBQgSUIAgAQUIElCAIAEFCBJQgCABBQgSUIAgAZ02+3OeYAUUSUAnza2gwJWATpl58+oLYNUEdMK8m+uXKd4MUDgBDc2roICABucVUEBAg/MKKCCgwXkFFBDQ4LwCCghocF4BBQR00ryGMQFdAjplZgPpgQ4BnTS3WzmBKwGdNrt+AhcCChAkoABBAgoQJKAAQQIKECSgAEECChAkoABBAgoQJKAAQQIKECSgAEECChAkoABBAgoQJKAAQQIKECSgAEECChAkoABBAgoQJKAAQQIKECSgAEECChAkoABBAgoQJKAAQQIKECSgAEECChAkoABBAgoQJKAAQQIKECSgAEECChAkoABBAgoQJKAAQQIKECSgAEECChAkoABBAgoQJKAAQQIKECSgAEECChAkoABBAgoQJKAAQQIKECSgAEECChAkoABB5QV0u2l9nj6zgALpFBbQ3ebW1IYKKJBOUQH9/WXz2sefk15BQIF0Sgro3++3xfzz7fDnT7+mvISAAumUFND9XS6bpL5MeQkBBdIpKaDb+wP2w07opNOgAgqkU1BAD7ubX+++uZt2DC+gQDoFBfSwt3l/uL6fdhlJQIF0BBQgqKCAOoQH8lJQQF1EAvJSUkD7hzHd75UOEFAgnZIC2juQftqtSAIKpFNSQE/FvPXhx6RXEFAgnaIC2n0Uk4eJAEsrLKD/PM4OyEZ5AZ3g7oB/I6BAOgIKEFR1QO8JKJCOgAIECShAkIACBAkoQFBBAe25D2nyx8oJKJCOgAIEFRTQvk81FlBgOSUF9LgPOunpdXcEFEinqIA2BZ34+KVXBBRIp6yANkfxkz7C4zUBBdIpLKD/dpueT5Ybb9zium8eGKO0gB4O4h/ZBR21uJ48AoxSWkAf3AUds7iXcCooMKi4gD5mxOJ2simgwBABHZpEQYEBAjo0iYACAwR0aBIBBQYI6NAkAgoMENChSQQUGCCgQ5MIKDBAQO8nMYwJGEVAe6YxkB4YQ0D7JnIrJzCCgPZOpZ/A+wQUIEhAAYIEFCBIQAGCBBQgSEABggQUIEhAAYIEFCBIQAGCBBQgSEABggQUIEhAAYIEFCBIQAGCBBQgSEABggQUIEhAAYIEFCBIQAGCBBQgSECZ2eZs6fcBjxNQ5rXZKCjVEFBmdQmnglIBAWVOnWxaF5RPQJlTdwVYGRRPQJmTgFIVAWVOAkpVBLR3KpeJn0RAqYqA9k2koM8ioFRFQHumMdDmaQSUqgjo/SQG2jyPf12qIqBDk9jGk7N/T00EdGgS23h6zjBTEQEdmsRG/gT6ST0EdGgSWzkwQECHJhFQYICADk0ioMAAAb2fxEAbYBQB7ZnGQBtgDAHtm8iFYmAEAe2dSj+B9wkoa+R3JEkIKCvkLA1pCCjr4zohiQgoq2OkGqkIKKvjXglSEVBWR0BJRUBZHQElFQFldQSUVASU1RFQUhFQVkdASUVAWR3DmEhFQFkfA+lJREBZIbdykoaAskb6SRICChAkoABBAgoQJKAAQQIKECSgAEECCtAxZYybgAJcTbrLQkBZIwPpecO0+3wFlJnl0C63cvKGiU+aEVDmlUO7PEyEt0x81qGAMqsc2uVxdrxJQIfYXhaWRbs8UJk3CegQm8vCsmhXFm+CPAnoEJvLwrJoVxZvgjwJ6BCby8KyaFcWb4I8CegQm8vCsmhXFm+CPBnGNMTmsrAs2pXFmyBTqxhI/+dbM5Dw06+p89lcFpZFu7IYCkCuar2V8xDNtpjbdgk/T3wN28vC8mhXDoNRyValDxO5BnR7/R0xcSfU9rK0PNqVw+1Q1KDEgO4PP/dfz9+ZuA9qg1lcHu3K4T1QgQIDevjvx5/n7+02H35MeQ1bzPK0i3oUGND9ef+z8ff7tF1Qmy2QToEBvdnr3E07CyqgQDoFBnTbDej+ejjf81Z6POvdAetTZEA7zRRQYDEFBnQ3eg/0noAC6RQY0P1m83L55tY5UCClygfSd6+8d6/IjyGgwKCKb+XcHIfQX6p5SOm0I3gBBQZV+zCRc0CPjtk83tH58u58XQIKDJj4tIaCAvqv09BjQHeT+ymgwJAVPFD595f2cvzUhzEJKDBkBQF9gIACAwR0iIACAwR0iIACAwR0iIACAwR0yNi7C9w3D6tU9TCmh428uUBBYaWqHUifwrh/khw+tQdYRK23cqYw6qRGDp8bCSyk0oeJpDDxrLCCAgMEdGgSAQUGCOjQJAIKDBDQoUkEFBggoEOTCCgwQECHJhFQYICA3k9iGBMwioD2TGMgPTCGgPZN9NZIWrd4Ah0C2jvVcD8VFGgI6JSZHdsDHQI6YV5Xl4AuAQ3Nq6DAswPa/Sj3juOHEi9CQIF0BDQ0r4ACAhqcV0CBuc6B7jebDz+OXzVJfUn9V44noEA6swT095fN18sf9m1Ll1BSQI05hdzNEtDt5vObf5pXQcOYjNqH7M0R0MNhe/eofV/oOdB5B9IbtQ/5E9BJc8+2U2jUPhRgjoD+/f7qEP7Tr9R/6VgPLu58B9WuWEEBZjoH2tkF3W0KPQc6KwGFAsx1Fb49bG++LPQq/KwEFAowzzjQ3c04+kLHgc5KQKEAMz1MZH/N54L7nwW1SEChALM9jel8V+dyF+CPimmRgEIBPM4uT4YxQQEENFMG0kP+BDRXbuWE7D3/cXYv9w+1K/VOpHnpJ+ROQAGCBJQ1sntPEs6BskJOMJOGgLI+hjiQiICyOgbZksq8AW0eJbLcs+z+2V5ouM2LVGYK6J9vTTebfi57M7zNBQElndk+lfPjz+bBykvfD29zQUBJZ7bngR72Ow//ObRzt+l8ROfcStpcXCd+FgEllVkCen4I/TmdpX4q58yMtHkaASWVmT4TqTkDevjP8fRnuR8qNycjbZ5HQEllxk/lPPznWE4BHcFImyfyj0sqMwZ0f/40OQEdwU7SM9m9J5EZA9p+Nueu3I81no+APpUTzKQx0znQr9dToIecuoj0LgF9Lv0kiZk+F/6Qzt3pJqRmSNNyH8tZzBYjoFCA2caBnj/P+HAcv+TNnMW0SEChAHN+LnwTzt2yH8xZTIsEFAow4+fCny/BL3f8/q+gFhlpAwXwOLtMGWkD+RPQXBlpA9mbLaDNByM11+IXfRxoSQE10gayN9vzQM9PAt0tOYiprIACuZsnoOfP5TwEdLtZtKACCqQzS0CbJyl/+vX7S3Mn0tYDlYFKzBLQ/fFBoKeAeqAyUIuZbuU83cR5DOj56aDLEFAgnfkeJtIG9LAL6nF2QA1mfJxdG1DPAwXqIKAAQQscwm+dAwWqMNNFpGaf8xxQD1QGKjHXMKbPbUCbMaEeqAzUYJ47kY6j548B3XmgMlCLWW/lPPFAZaAOMz1MpDlwz6CfAgokNNvj7M4JXTSfAgqk5IHKAEHzDGNaeL/zSkCBdGYaSL/sc+ivBBRIZ6ZbOZcbOn9LQIF07IECBM11J9KinwZ/JaBAOvNchf/9JZN90AcX1wdlAh0znQO9Verj7HxUO9AloFNm3rz6Alg1AZ0w7+b6ZYo3AxTOnUiheRUUENDgvAIKCGhwXgEF5g3ovjn/uew9SQIKpPP8gP7+choCenkk6JJj6gUUSOfZAW2yeQrotr0Gf/pszmUIKJDOkwN63O08jlr6/eV0+L4t9jORDGMCbj05oPvLIfuusyda6KdyGkgP3HhyQLftAfthV/T81Z9vm6+p/87R3MoJpPPcgB5ieT5ev351iOpDV+Iv51IjL+NhIkA6Tw1oc97zDcHzoLtXLzO1odoHpFNUQHteb+Jt9QIKpPP0Q/jT+c6/36+p20YfJnK9pH959aklFlAgnecG9NLN/fVoO/4JSfu7XDZJnXRJX0CBdJ58FX532mU8HHu3w+eb5gWvwvfsuk6tsYAC6Tw5oNdHgZ47t48/DvSQ3vvy7qYdwwsokM6zb+VsC3qK5ukqUHAc/eGl7ufcD+W47+pV7O8GuPf8h4kcRx6ddx2POY3ehySgQF7mfR7oQ3chOYQH8lLSA5VdRAKyUlJA+4cxTdqlFVAgnZIC2juQftolfQHlyClxkigpoPefjzz56cy2mOXl0C4XFUmjqIB2H8V0M7p0NBvM4nJo1+UvV1AeU1hA/y37ODselkO7On+1HwgeUl5AH2J7WVgW7fLpVqTy9Kcx9Qney5nizdlclpVFu7J4E1RBQJlTFu3K4k1QhScfwvc/Ujn3gC5/laNaWbQrizdBFeZ4mMhynyF3Z9Ti5nCduFZZtCuLN0EVnn4R6VDQiWM1n2nM4uZwnbhaWbQrizdBFZ5/Ff5wFB/8BLknGLG4WVwnrlYW7bKKSWWGYUy7+BPskhsV0EmTM0ke7XKQQSIzBLTzkfCLE9Cl5dEup7lJY46B9Bntggro4vJoVw7vgQq4E2loEhvYE2gX9RDQoUls5cAAAR2aRECBAQJ6P0kO14mBAswR0D/fbobSbxe8KG8gPZDOIgHN/F74PK4TA4uYsvnPH9A/33IPqOvEsF6TdqCeHNBd7+Ps8j6EB9Zr2im8Jwe0/4Ggyz2eSUCBARMvIj/7EH6fVT8FNANOkJCxicMYF7iItCSb7eJcoiNnAjrEVrs0g8TIWoYBzYiNdmFuUyBvAjrERrswN8qSt6wC2lyEf7m/FJ/7OFCeRkDJm4AOvjnb7LIElLxlNYxJQLkloGQuq4H0ubHNLkxAyV1Gt3L+/Z7PxyEd2WYXJqBkb8pA5acfwn82DpQrw5ioioAyKwPpqcnTD+EP7RRQrtzKSUWefBFp2/cwJlfhV00/qceTA/r7i4ACtXr2MKbeggooUANPYwIIElCAIHciAXTkM5A+OwIKDMroVs7sCCgwxMNEBggoMCCrx9llR0CXZyA9GcvqgcrZsdkuzq2c5ExAh9hql+ZhImRNQIfYaBfmcXbkTUCH2GgX5oHK5E1Ah5S0zVZ5qlBAyZuADilom63zYouAkjfDmIaMvT1r+XJVerFFQMmcgfQDRt6dtXxBa73YIqDkzq2cbxv3OyWDfb9aQ1PrclERDxN506izwjns+9Uamjz+dSERAR2aZLltPI938QRZ7N9DIgI6NImAppfDGWZIRECHJhHQJ9BP6iGgQ5MIKDBAQIcmWTKgLrZA9gR0aJJFL8O72AK5E9D7SfLY93OxBbInoD3T5LHvp5+QOwHtm0i7gBEEtHcq/QTeJ6Cskd+RJCGgrJCzNKQhoKxPLtcJKZ6Asjq5jFSjfALK6mRyrwQVEFBWR0BJRUBZHQElFQFldQSUVASU1RFQUhFQVkdASUVAWR3DmEhFQFkfA+lJREBZIbdykoaAskb6SRICChAkoABBAgoQJKAAQQIKECSgAEECChAkoABBAgoQJKDMzE1A1ENAmZfb0KmIgDIrD0KiJgLKnDyKk6oUGtA/35pjwE+/ps5no12Yh8FTlYICeohmW8xtexrt88TXsM0uTECpSpEB3V4vREzcCbXNLkxAqUqJAd0fwvn1/J2J+6DjFtdl4qcRUKpSYEAP//348/y93ebDjymvMWpxDbR5HgGlKgUGdH/e/2z8/T5tF3TM4hpo80QCSlUKDOjNXudu2lnQEYtroM0z+delKgUGdNsN6P56OD/GqIBOmpxp7N9TkyID2mmmgJbGGWYqUmBAd/ZAi6af1KPAgO43m5fLN7dPOAc6aXJgvQoMaPfKe/eK/BgCCqRTVkA3xyH0l2oeUjp4BL/p8e5fI6DASOUF9OiYzeMdnS8DcwQDaqANMEpBAf3XaegxoLt3+tnDQHognbICevL7S3s5furDmNzKCSRUYkAf4GEiQDoCChBUbECnPofpRECBdAQUIEhAAYIEFCBIQAGCBBQgSECzZTQq5E5Ac+V+KMiegGbKHfmQv2IDGlNMjDwTCgogoHnyVFIogIDmSUChAAKaJwGFAghongQUCiCgeRJQKICA5klAoQACmifDmKAAApopA+khfwKaK7dyQvYENFv6CUuYsuUJKMDVpGM/AQW4mHb1QUABWhPHvwgoQGviCGwBBWgJ6BABBQYI6BABBQYI6BABBQYI6BABBQYI6BABBQYYxjREQIEhBtIPEFBgkFs53yagwDAPE3mTgALpCChAkIACBAkoQJCAAgQJKECQgAIECShAkIBOm90HZQIXAjppbgUFrgR0ysxTnjIAVE9AJ8w76TlXQPUENDSvggICGpxXQAEBDc4roICABucVUEBAg/MKKCCgwXkFFBDQSfMaxgR0CeiUmQ2kBzoEdNLcbuUErgR02uz6CVwIKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBJUX0O2m9Xn6zAIKpFNYQHebW1MbKqBAOkUF9PeXzWsff056BQEF0ikpoH+/3xbzz7fDnz/9mvISAgqkU1JA93e5bJL6MuUlBBRIp6SAbu8P2A87oZNOgwro8tqzL0u/D3hcQQE97G5+vfvmbtoxvM12cdfz10u/E3hYQQE97G3eH67vp11GstUu7RJOBaUCAsqcOtm0LihfQQF1CF+B7gqwMiheQQF1EakCAkpVSgpo/zCm+73SAbbZhQkoVSkpoL0D6afdimSbXZiAUpWSAnoq5q0PPya9gm12YQJKVYoKaPdRTB4mUiQBpSqFBfSfx9mVzTAmqlJeQB9io12agfTURECZl1s5qYiAMjP9pB4CChAkoABBVQf0btSoA0cgIQEFCCoooD33IU3+WDkBBdIR0CnsxgIdBQW071ON5w2oEwFAV0kBPe6DTnp63Z2HFtdNNPXwm5AkigpoU9CJj1965ZHF7WxuNrzCOZYgjbIC2hzFT/oIj9ceC2ia12FxjiVIpLCA/tttej5ZbjwBxbEE6ZQW0MNB/CO7oAKKNUk6pQX0wV1QAcWaJJ3iAvoYAcWaJB0BDc1rsyuZNUkqAjphXpce6iCgpCKgU2Y2+KUKAkoqxQZ0FxpS71ZOHEuQjoBOm10/a+BYgkQElBVyLEEaAsoa6SdJCChAkIACBAkoQJCAAgQJKEBQsQGNEVAgHQEFCBJQgCABBQgSUIAgAQUIElCAIAEFCBJQgCABBQgSUIAgAe2dytN2gfcJaN9ECgqMIKA90/jEMWAMAb2fxGfeAqMI6NAkCgoMENChSQQUGCCgQ5MIKDBAQIcmEVBggIAOTSKgwAABHZpEQIEBAno/iWFMwCgC2jONgfTAGALaN5FbOYERBLR3Kv0E3iegAEECChAkoABBAgoQJKAAQQIKECSgAEECChAkoABBAgoQJKAAQQIKECSgAEECChAkoABBqwsoQDrJG5X6BVNa+h8bqEvyRqV+wWKs4PyARayBRcxZue/8UQWvtLEsYg0sYs7KfeePKniljWURa2ARc1buO39UwSttLItYA4uYs3Lf+aMKXmljWcQaWMSclfvOH1XwShvLItbAIuas3Hf+qIJX2lgWsQYWMWflvvNHFbzSxrKINbCIOSv3nT+q4JU2lkWsgUXMWbnv/FEFr7SxLGINLGLOyn3njyp4pY1lEWtgEXNW7jt/VMErbSyLWAOLmLNy3/mjCl5pY1nEGljEnJX7zh9V8EobyyLWwCLmrNx3/qiCV9pYFrEGFjFn5b7zRxW80sayiDWwiDkr950/quCVNpZFrIFFzFm57/xRBa+0sSxiDSxizsp9548qeKWNZRFrYBFzVu47B1iYgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShA0PoC+ufbx5/t13+/by5elnxT6fz+0izMhx/X7/z51nzn63JvKbXXi1jtWuyus0rX4nWBCl2LqwvoYT1dA3r6qSxupb3tukCXhdyev9FtasnuF7G6tXhtyWVxaluL94tY6FpcXUC3nbb822+KXGlv6v4Qnpdy9/obhetZxNrWYs++WG1rsWcRC12LKwvoccVdfwZ3Ra2r9zXb2fGgqNlh+dx80Rwoffp1+vH8vOybS+N+Eatdi83vimbdVbwWr4tY6lpcV0BP512uAd1W8gv9rPlxvO6yHBdt225xh0Wv4fCvZxFrW4uHX/LnNXVY2NNXta3FnkUsdS2uKaDHo7//7JwDPazH86+/Ouyv+yfNrvbLcZHbpd1VcQnifhGrW4uHSLaLeN4rq24t3i9isWtxTQFtjhteuheRDj+YVRwQtbpb1/b4g9n5ST18WeRP6K37RaxuLXacl7C6tdhxXsRi1+K6Avr59ir8YXfm5bhbWsNh0Svb40J1enMz/KAKp0WseC22e54Vr8XLznWpa3FNAT3q/gTuNh/+63zdr65f68cfzGaROmfmr+edKnFexHrX4q69Hl3vWrwsYrFrcdUBbQfXFbfW3nfeadneDCWsadO77pfVuRa3nfE8la7F7iIWuxbXHNDmKsTpp3G/qeLc/FV7qqy7uVW06TXaRax0Ld4GtMq12F3EYtfimgN6HULRbI01nVq6LM7tplfkQLs3XBaxzrXYjjTv+TVYy1q8WcRi1+KaA9pVz8/lv+Mv8fMyVrrv0l3ErqrW4mncXXNputq1eF3ErqLWooCe1DG+7qQTl1o3vf5+VrUWG+f9slrXYqOz69kqai0K6Mm+pJU2bNs5C1/p9dvtGxcaKlqLJ6eWVLoWT+5zWdRaFNCTolbakObM0vWYaF/jCMLbReyqZi22TgtU5Vps3a+zotaigJ4U+iiDO83d/l9v/ljdPSyvFrGrlrV4cWpJjWvx4j6XRa3FNQf09gaPKo6MXj9sorq7qO8Xsbq1eH+3anVr8X4Ri12Law5oZ7zEvo7HhN2fka/tOT73i1jdWuwsUPtlbWvxfhGLXYtrDmhzLu309a64W3B7tc8n6qjtSZL3i1jlWrws0Gm/rMa1eLuIxa7FNQf0/HjQkwoOjG4e6t25j/qsiosPPYtY21q8WaBzLitbiz2LWOpaXHVAO2utnLPWb+t+TkL3PsCTkn6tv6l3Eetai/+6H1ty2d2sai3+61vEQtfiugPa7tHUcFh0+3FBnR/Emj7P8Y1FrGktHp1q0r3eXtNaPLpfxCLX4uoCCpCKgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAUqK/3z/+XPo9gICSid3m4NP/fd+cy/jnf7r/98+3bjD3zbSbDz9GzH37Mj3+fNt8+vXom2etBJQc/P1+bOLnv20Ct5uvt//7GtDfXzZnbfnenvv2ZfoIKA8QUHJw2qfcfG0TePhzt3zbzTWg135eCvrm3K9epo+A8gABJQe7zW3Gbsp33MG8BPQQ081Lcw60OWx/GZ5bQHkuASUHu7aFZ93ynXY524Aegvfhx+ki0r7t5ptzCyjPJaDk4O2AHgK32fzn9RzoIaeH4B0DejnlKaAsREBZ3LGRx5OY5yaez2keq7g7HbH3BfS9ubsv054ofbn9Wz8LKI8QUBb3TkA/31yF7xzCvzd392Uul57a+Xbny1D/K6DECSiLGwzoUXcY0+50vX1aQDuX7i9H/Sf/IaDECSg5aM9ivjGMqRvQ86DP7sD6N+duv2jmOY673x8P20/Vbbq57wyGgskElBxMCOh177G7Uzoc0Mv1+vMZgOvAp2bXVECJElByMCmgl5Hzx53JobnbL7bX+z7352P99huvx5DCBAJKDiYGtD2OP0/0XkC7V9oPu5zHS++3V/UhREDJwfSAHv68bR8oMiKgXZ9+dYpqGBMPEFByEApoU9Cvg3Ofv+jePn8MaGe3U0B5gICSg1hAj4fjQ3NfA3obSXugpCGg5Eue3bkAAAGaSURBVGB8QLfHKacFtHPK88Q5UNIQUHIwPqC7zenWpObPu5GH8M0lp5ub5TvfcBWeBwgoORgf0MMeY3srZ3Nq82Vw7v21sO38p29dstmOqIcIASUHPQnsZq17DnTb3sq537x+nN3d3O0Xx8vwL+e5myv3zTfOJwEElDgBJQevE3i6bt57L3x3TNJ5NPybc1++2Hcuwrf3Jp38t4tIxAkoOXidwNNA+ctR/M1V+GtBO3cT9c99fZnrSKbzi54L6mlMPEJAycHrBJ7S97n936+GMZ1r+PX9uTsvs+ve/Hn9X4Yx8QABpUQ+F54sCCglElCyIKCUSEDJgoACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKECQgAIECShAkIACBAkoQJCAAgQJKEDQ/wNjkD9WcsrN3AAAAABJRU5ErkJggg==)
qqnorm(fit1$resid)
qqline(fit1$resid)
![](data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAABUAAAAPACAMAAADDuCPrAAAAclBMVEUAAAAAADoAAGYAOjoAOpAAZpAAZrY6AAA6OgA6Ojo6kJA6kLY6kNtmAABmAGZmZgBmZjpmtrZmtv+QOgCQOmaQkNuQ2/+2ZgC225C2/7a2/9u2///bkDrbtmbb25Db/7bb////tmb/25D//7b//9v///8dkKWdAAAACXBIWXMAAB2HAAAdhwGP5fFlAAAgAElEQVR4nO3d6YLayMKYYWYysbPYST5PMiQ5JBns4f5vMSAk0AqotFbV8/w4bveo1ZJB79FSEocLAEEOWy8AQKwEFCCQgAIEElCAQAIKEEhAAQIJKEAgAQUIJKAAgQQUIJCAAgQSUIBAAgoQSEABAgkoQCABBQgkoACBBBQgkIACBBJQgEACChBIQAECCShAIAEFCCSgAIEEFCCQgAIEElCAQAIKEEhAAQIJKEAgAQUIJKAAgQQUIJCAAgQSUIBAAgoQSEABAgkoQCABBQgkoACBBBQgkIACBBJQgEACChBIQAECCShAIAEFCCSgAIEEFCCQgAIEElCAQAIKEEhAAQIJKEAgAQUIJKAAgQQUIJCAAgQSUIBAAgoQSEABAgkoQCABBQgkoACBBJSRfn49XH17/P3X9+tff/tr+V/8+hfdF+vqx9sZvZ20+E01Pxr/4f26nr+9m4JECCgj3fPzrMgeAvrPn/Xe/f6vF3P5ZNJ2QB8J/Whdr/9AApoLAWWkcv/tS/X3HQT01O7dcMA+mrQb0LKgH6xrEWgBzYWAMlJ1AFwd124f0GMnd8+8B03aE9D77/1gXY+vA05aBJSRqoD+8ff975sH9FRPXLV0/QX9cNLGbzrVphJQmgSUkR6XYMpKbB3QcnkezTo1d5BDJm3+pvtPFWdLBZQmAWWkR0DLCzAbB7S8KFSL4Lmxgxwyaes3HR8/JqA0CSgjPQJaHvu2o1JmqnZ9+3RvSvlHkZjb9Kdaztpjoy7Pyz2Pb/bX6/6TjcPwY23GQZO2ftN5MKDtdT0fnt4PpyJ+AspIz4D2XVh5/tfn9+7lPFZduQf0cTn8uv/3uGjzbFP9avnLcwXHQ+e7PaEcN+mHAe2uq4DmRkAZqcjG7/+1bF//FZdm+Ypv/sdDtbNWBPQ/PKf64/8+Q1Tty/XNpzeg933A5kF43/dGTfpZQHuWUUBzI6CMdA/o//reF5Xz4dCtSntvsmcoUfsnavt2z6z2BvS+89rahTzWWxwwad850OKv79ZVQHMjoIx0D+i/TlV5alEpD8WfKbl/uxbQ2zfKgH6p/ZfbfOoXdI6Nmb8ahNlz8rSc7cDZ0k8mbf6m8zPi79fVRaS8CCgjlQG9H/x+a0SlfpKxtr9XH0lZTVU25lir5qmvU/fq/bhcNgvo8bmD+cG6CmheBJSRyoDe97uuBXlG5d6RKh7nR2FOtUpemkfN9Svx59r3r3muXeN/EdDzoSdY/QH9fNKeO5HaZxGG1lVA8yKgjFQF9J6KL7Wo1AtYXZ+5paReyUv5c2VOT7UfKWbc37It9kAb6ruaL9dVQPMioIz0CGh5dP2Myqmxo1kFthOpx/cvzQ71BfQx9umyaUCrUr5fVwHNi4Ay0iOgZUT+XzMqz6vcj7/WdzMvPQEtO9QKaP3a/XBAX19ar13N/zb2KvzTs5O/3q6rgOZFQBnpGdB7Z/7TZwHt21sbDmhrjNBwQO8Hz80C3pfrNtdGQF9P2vlu3x2bAkqTgDLSM6C1zL0/hB8V0Gqo6B9/vzsH2rgZ89f32nD2b5dWQF9PWvdpQB3CZ09AGakW0OfT3QcvrNxyNTqgtYvabwNavxnzVHx1/8XFhM2Avpy07oOADq2rgOZFQBmpFtBnoN4PYxoT0FqE3ga0jPi3x7S//bfevcoxk34QUMOYuBFQRqoH9HGs/X4g/YiA1gYF1YdADWTt/NjDrF/86f2wo08n/SCgBtJzI6CM1Aho1aHGM957b+UcHdCimm+HMV0el+sbx+wDz+z8cNJPAjqwro1RriRPQBmpEdDmh2R0PrLtx/O7Yw7h+z/5bcxnIg0l7LNJPwnowLo+v20vNAcCykjNgJZ7i+8fZzcmoPVxmP/uazX1iE/lPBwGDuI/m/SjgPava+cTT0iagDJSM6DN49cXD1QeNYzpOZf7jU7vPpCo+WHvLwv6yaSfBbR3XZ87uUMfDEpKBJSRWgHtPOd94CM9xg2krw3QfFyVef2BRLVTmrWr7WGTfhjQvnV9zH5wQUmJgJKgW8P6d0CnTAptAkqSzp+fgxwxKTQJKEAgAQUIJKAAgQQUIJCAAgQSUIBAAgoQSEABAgkoQCABBQgkoACBBBQgkIACBBJQgEACChBIQAECCShAIAEFCCSgAIEEFCCQgAIEElCAQAIKEEhAAQIJKEAgAQUIJKAAgQQUIJCAAgQSUIBAAgoQSEABAgkoQCABBQgkoACBBBQgkIACBBJQgEACChBIQAECCShAIAEFCCSgAIEEFCCQgAIEElCAQAIKEEhAAQIJKEAgAQUIJKAAgQQUIJCAAgQSUIBAAgoQSEABAgkoQCABBQgkoACBBBQgkIACBNp3QA8As5k/UbPPcUZb/2sDaZm9UXPPcE4L/B8GkJNGNgUU4FOtvU4BBfhQ+6BdQAE+0j3nKaAAH+i7ZiSgAG/1X3IXUIB3BkYsCSjAa4MDPgUU4JUXA+YFFGDYy/uNBBRg0OvbNQUUYMC7u90FFKDX+6eFCChAj08etiSgADVjnlUnoABPox71KaAAD0U4Py6ogAJUbtks4ymgHQIKvFD2s/r6g8nnXoC5ZzgnAQVeqB+6C2ibgALDxn7gkYACFMZ/4JGAAty0L70LaJuAAr26l94FtE1AgR7VzufIz3wXUCB3tWP37hevf3D2RZl7hnMSUKCtcerTrZzDBBRoardyRD8FFMjZx60c+OkZF+U+x7lnOCcBBZ6m5VNAgXxNzKeAArmanE8BBfI09ei9nMkci9KY49wznJOAApeZ8imgQIbmyaeAAtmZK58CCmRmpqP3cl4zzeg5x7lnOCcBhbzNmU8BBXIyaz4FFMjHzPkUUCAX8x69l7Ocd34CCuzRAvkUUCALS+RTQIEMLJNPAQWSt8jReznn2ec49wznJKCQm+XyKaBA2hbMp4ACKVs0nwIKpGvJo/fyF8w+x7lnOCcBhWwsnk8BBRK1fD4FFEjSGvkUUCBBKxy9l79n9jnOPcM5CSikb618CiiQmtXyKaBAWlbMp4ACKVnv6L38dbPPce4ZzklAIWEr51NAgWSsnU8BBRKxfj4FFEjC6kfv5W+dfY5zz3BOAgop2iafAgrEb6N8CigQu83yKaBA3LY6ei9/+exznHuGcxJQSMqm+RRQIGLb5lNAgWhtnU8BBSK18dF7uQyzz3HuGc5p839vYBZ7yKeAAjHaRT4FFIjPTvIpoEBs9nH0XhBQICY7yqeAAlHZUz4FFIjIvvIpoEA0dnX0XhBQIAr7y6eAAnHYYT4FFIjBLvMpoMD+7fHovSCgwL7tNp8CCuzcfvMpoMCu7TmfMQb0eKh8Gf/De34pgLYdH70XIgvo6dA0tqG7fi2App3nM7KA/vx6aPv9X6PmsPNXA3jafT7jCug/fzaL+ev79e9//D1mFrt/PYC7vR+9F2IK6LmTy1tSf4yZxf5fEOASST7jCuixe8B+3QkddRo0hpcEiCOfUQX0urv5rfPN07hj+DheFMhbLPmMKqDXvc3u4fp53GWkWF4WyFckR+8FAQV2JKZ8RhVQh/CQvKjyGVVAXUSCxEWWz7gC2j+MqbtX+kJkrw7kJK6j90JMAe0dSD/uVqTYXh7IRoT5jCug92I2/fbXqDnE9wJBHmLMZ2QBrT+KycNEIB1x5jO6gF48zg6SE+XReyG+gI7QOeCP9VWChMW8YQoosKWoN8ukA9oV8SsFKYo6nwIKbCf6o0IBBbYRfT4FFNhI/PkUUGATKeQzqoD23Id0GHsvZwovGcQvgaP3goACK0sln1EFtO9TjQUUopNMPuMKaLEPOurpdR3JvG4Qq4TyGVlAbwUd+filloReOYhROkfvhbgCejuKH/URHm0pvXQQncTyGV1AL6dDzyfLfS6tFw/iEp7PvT7LIraAXg/ip+yC7u8FgFxMCOBunwYUW0An7oLu7t8fMjElf4+f3F1BowvoNHv754c8TNp7rP3o3rZgAQWWNu3gu/6zO9uEBRRY1tRzlwK6Fzv714f0Tb/2I6B7sbN/fUjdHJfOBXQvdvavD4mbZeSRgO7Fzv71IWkzDdwU0L3Y2b8+JGy2ge+GMe3F3v75IVVz3jdkIP1O7O2fHxI1722XbuXch939+0OKZk/dTvspoMDMdpm6ZQgoMKeM8imgwKxyyqeAAjPKK58CCswmq6P3goACs8gvnwIKzCPDfAooMIcs8ymgwHQ5Hr0XBBSYJtt8CigwUb75FFBgkpzzKaDABBkfvRcEFAiUez4FFAiVfT4FFAgjnxcBBUI4ei8IKDCWfJYEFBhJPisCCowin08CCozg6L1OQIGPyWeTgAKfks8WAQU+M5DPzz+zfa+f7h5OQIFPDKXvcPi0i59PGQ0BBd4b7N7j22+7+PmU8RBQ4K3B3cba999sXZ9PGREBBd54cdRd/w+vN6/Pp4yIgAIvvTxrKaBzz3HuGc4pmdcNVvLmoo+Azj3HuWc4p2ReN1jHuwvnAjr3HOee4ZySed1gDYeaoSn6v54yZUQEFOhXD+fwVfj+r6dMGREBBfq0djwNY+ojoECPIp4f7DUaSD/zHOee4ZzSeeFgSeW+5yeH3W7lnHeOc89wTgm9crCYR+Q+Om/pYSKzznHuGc4ppZcOltG4dFT/9jaLs2sCCtQNXTmy8fQQUOCpeYQtoG8IKFBpn6FMcujRnAQUuOu5wJPi0KM5CShQ6L0+nuDQozkJKHAZfuanfr4ioJA7jQwmoJA5R+nhBBTy9mingo4noJC1gfuO+IiAQs4+eWQdg9YN6M+vh8Mff8/9G0fwBoG64RuP+MRKAf31/dbNWz8Ph9/+mvtXfs4bBJ7al45sH2OtE9Dz4fD7vy7//Hl/va5fbsUbBCqHT5+ZzKBVAnrb87zud17/uLbzdDh8m/t3fswbBO6e197r39tscSK1SkCvzfxy/+OWzmPxl214g0DhcewuoFOsEdDrofvtDOj1j+L053nDY3hvELg0rh25CD/FGgH99f3w4/5HUU4BhU01Lx154NIEKwb0fD+QF1DYUuemTbdyhlsxoMdD8cfltOFQUG8RctcTSv0MttI50G/PU6DXnLqIBNvQyXmtchX+eEvn6X4T0m1I04+5f+fHvHnImf3Mua02DvTmR3Ecv+XNnN495Es+57fOnUin4qW7hfO06Y1IAkq+5HMBK90Lfz48LsFvd/x+EVCyJZ+L8Dg7SJ+j94UIKKSuzKfRSvNbLaC/vhcPFDlt+jhQASVDz0/sUNC5rfY80PJJoKctBzEJKFlohLL6wh2bS1gnoPd+3gJ6PGxaUG8dktK7T1nf1Tx0+mkrmNMqAb09SfmPv39+vd2JdPRAZZhJ71F5fVfzUPtLfYrVljB1qwT0XDwI9B5QD1SGmfQelTd2Nft3O20Gs1npVs77TZxFQMung27DO4d09Oex/qDPvm9fbAYzWu9hIlVAr7ugHmcH0/U3sXbtSEAXt+Lj7KqAeh4ozOFFQB+XkN5MzEQCCpEaDujz0vubiZlog0P4o3OgMIOhgPZeejeMaRErXUS67XOWAfVAZZjFi4D2TmEg/fzWGsb0pQrobUzo9JH095H54/dkvXVIR19AD0NDl9zKuYh17kQqRs8XAT2FP1D5Gs3qJ4/VO2Hsvqz3DunoSWVr7Hz/HUqrLWAGVr2V8y70EtIzoMfnzEbG2JuHhLRTWdbRruZ6VnqYyO3AfWI/nwG9PZ35W/mdkfug3lKkpJHKZzL1czWrPc6uTOiEEUxVQK9/PuZyun/S5+cL5z1FSp6pVMxNRPRA5Sqg59rN9Ncsj9oF9RYjTfK5jQgD2tjrPI07C+pNRorkcysRBvRYD+jL25oOPZZaOtiKN/Z2lg1o4+r71Ovwz4DWfl5AyZy39ZYiDOjp4z3QnoXzTiMt8rmpCAN6rt/KNPLGeu81otV3ECWfG4vwHGj9yvt55OPtvduIVc9pKEfvm4sroDffntW8jS0dtzfr7UakHqmsf+H9vLX4Avo8CVDc0TnuwSTecMSp1sra7UabLQ6liAJ6qTW0COhpdD8FlEi1n7wkn/uw/EWkH91LSROfSP/za3U5fvSDRb3piFPryZ76uRMxBnTKwnnXEaVOQLdbFGoEFCLQuva+4ZJQF9c50Mm884jTwKe9szEBhQgMfNo7GxNQ2KH2oPnD85mf3sQ7ss7nwjeeeuxjjeGN7m1Hj7/aAd2TTQLqIhK80r3tyCcd7dP6Aa1/IsfqvPmIQPe2o/qjGbdaKnosHNDToY9DeHih84HvwrlXCwe0/3l2ox6gNCvvQiLgvs1oLH0If95VPwWUGHTv29xuWXhpg4tIW/JOJALu24yGgMLetC69b7kovGYgPeyN+zajIaCwN2U1nf7cv5UC+vNr/SqScaDwyv2OTQHdv3UC2hoOKqBQ14mlYfORWCWg7bFMAgo13VzqZyRWCehxy2Y2eEOyP49OVl8oZzTWCOg/f+6lnwLK/jR2PC8uvUdlnXGgW9581OCNye74vLiICShsy31HEVvpEF5AYYDPi4vYKheRThs+wK7J25Pdcd9RxFYJ6H6O4b0/2R2fFxexdQbS//y6k31Q70/WNzhMvvpb9T1v0OisdBHJQHqy1Rkm3/rGobxv8+L9GR8BhUX1DZPvfKPqp/dnZAQUltQeJt/7DXduxsrj7GBJPR8Q1/yGfsZMQGFJbwIqm3ETUFjS+4CuvkjMZ4uA/vovzoGSixcBlc/4rRTQ5hNBXUQiG8MB1c8ErBPQo6vwZGoooIf78M9tFoq5rBLQ5iciHQ4/5v6dH/OGZWX9w5gO+pmGtZ5I/+X2SJEv96+366d3LKsbGDd/MG4+BSs9zu52J/z5fuh+2vLx9N6xrK4zzNO4z3Ss+EDln19/++tS5NQhPDlp5tK4+ZSsFNAfzz/KQ/lteM+yLeFMy4oBrR5M//Prdo+2895lcS/2L+UzNSsG9HK873r+/GoYE+l6cYQun8lZ8zORyk/2OBsHSro6F93r/8XbLzUrfSZSkczzobiKtOUnJHkHs6zOsM/af/DuS89aA+lv6bx/sMf54CISyerceFR9KZ9JWu9Wzm/POzoNYyJV/QGVz1StE9B//iyqWT6afrsdUAFlYX0Blc90rfU4u3Nx+rPaFd2MNzKza1x07wbU0XvKPFAZJmkOW+p5+JI3XcIEFKZoDVtqBVQ+EyegMEF72NKheeVIPxPnY41hgr5D9upv8pk+AYUJBi4a2fvMhIDCBD3DluQzI2ufAy0frrwVb2rmNTxw3lstB+tfRNryVngBZWZDA+e90fKwfkA9kZ6EDBzCb7Q0rG2DYUyexkQ6Ok9fks+sbBJQF5FIRnMgvaP3zAgoTFIbtCSf2Vk/oL++O4QnJY1+br0wrGv1gN4GhXqgMumRzxxtMZD+/mS7TXiLswxH73naIqAbPhDUe5wlyGeu1g/odleQLgLKIuQzWx5nB9PIZ8YEFKZw9J41AYVw8pm5LQL6z383kJ4kyGfu1gnoubyA9KP621aPE/F2Zz7yyQoBvX0o/MOX4u+bXYn3hmcujt5ZI6DnQ8Nvf503vBXJO555yCc3iwf0XBv7WQ0I3e5WJO95ZiGfFJYO6M+vjaHz98P57W5F8q5nBvJJaemAHg/NZy/dCrrhvUje90zm6J2HhQN63QFt5vK2R+ojPYiXfFKzcEBP7Voey0vxG/HWZxr5pG7ZgHZGLF13QP/NA5WJlXzStGxAO0+fPx1+/98bDgMVUCZw9E7b4gFtHa4fv1yP4gWU+MgnXWsH9HIRUGIkn/RY/Bxo53xn3/dWYxsgiHzSa+WLSJeBvdK12AoI0Dx6r+5K3nCB2Iu1hzEV93a6E4mItGr5fK7DdovEXqw9kP62A+pTOYlIK5WPvykoq9zK+eX1N1blPc847T3N2l+9mVjlYSK1a0Yb3wrvPc8o3SP15sH8yovD7qzzOLsyocXj7DY8gPeWZ4y+E50CSt3yD1QuHh9Sf6Dy3L9wDG95PtZ7nUhAqVvjM5GOzY/02JK3PB8auMwuoNSt86Fy5V7otnufN97yfGRwmJKAUudz4aHtxShPAaVOQKHl1SB5w5ioE1BoeHOPkYH01Ago1Ly/R9OtnDzFF9DnNf2AK/re9LzyURf1k4fIAno6NI1tqHc9L8giI0UV0NaY/JuR94XaPhgkn4wWU0BvN9LXi1ncGTru4cy2EAY4KidATAE9d3J5S+qoD5m3idBLPgkSU0B7Pkxp7OPtbST0kU/CRBTQ6+5m91H2p3HH8DYTuuSTUKsF9HbC8re/LqfwD5S7zqF7uH5+dRmpc8nJhkKHtwXhVgpocb2nCOjIk5bNeQgoM/OmYIp1Anrv5y2gt2HwgQV1CM/s5JNJVgno7WL5H3///Hp7nN0x/DM9XERiXvLJRKsE9P5JxveA3g7iAz/WuH8Y06iZ2V54cPTOZKsE9FiErwzoNXqBF5J6B9KP2521wVAakU/nzxmyRkDLk5dlQK+7oKHH8OWp1LqRD7m3EXD3KoitYLoCyaA1AlpePq8C+vLC+RvHVj89TIQQL2vYCuZjUgWlI7KAXjzOjsle7022glmb1LuHtg0O4Y+h50BnYBPgzcF4O5g+BIkXVrqIdNvnLAM6duTRrGwB2Xt3LrMdTAHlhbWGMX2pAjr6AUqzsgVk7v2lIAFlhHXuRCpGzxcBvT1SfrsjeFtA3gaO3lvXjOr/QUB5adVbOe8mXUKayBaQsaGTn82r7gLKCCs9TKQYA799P20BGRs6eu9cda//JwHlpdUeZ1cmdNN82gLyNXjy891Vd8OYeCGiByrPwSaQpxdDl94G00B6hgkoyXs58rNziN4Jpls5GSSgpO7dwPn2151g6idDlg1oz9M/tj0RaiPIztiB8xfB5HMCSsreh9BVdiYQUNL1yX6kgDKBc6Ak66PDcMOUmEBASdSnZzENUyKcgJKkV0fvAxfZvTkYTUBJ0JuRn4YpMZPVAvrza/Ee3fBJTDc2kiy8rKFDduazUkDPe7gEfxHQLLzemXTRiBmtE9BTfRDTyA/SnJVNJnm1d1rne/ev69OuvnikZZWA3g7fyx3P2pdbsMUkrjnguPPNi4Ayq1UCeqx/isd1b/Tb3L/zY7aYtDXCWX3R/IaAMqM1AvrreyOZPpWTZRTt7JzjbH1DQJnROgFtnPac+Lnwk9hi0lXue/Y8nq7xDQFlRgJKGgYvEgkoy1npHOiXwb+tyxaTqPqZz/p3u98wjIkZrXUV/nkS9Oxz4ZlZ85J7/fs93zCQnvmsMw70WtDyKP72gLvtLsILaIraYz4vza/7PrPDvZvMY6Wr8H22OBNqm0lPq4TvA+red2YjoEStE8J3w5hgRgJKxPp2JN8MpIcZeZwd0Ro4Du+c43TSk6UIKLEaLGInl/rJQgSUOAkiOyCgxMgOJbuw/gOVN7p8dGejS4J8shOrBPSfP7e//n5ns0uBfLIXqwT0dBBQ5tLNp2tEbGWNgN52QLd7fkiDjSx2PaU0SonNrDSQfuMP43ywjcXto4HzsJoNnki/JZtY1Pr2Mmvf8uqytpUO4QWUyap8Ng/YPSCZDa10EWknp0BtYfE6tPv5uNO9Ps0mi0a+Vgno9Rh+u2coN9jCIlUPZvOUp4CyodUeqLyPg3hbWJwOfVeKBJTNrXQn0rWgxoES6NCpZu1rAWVD7kRi595cMxJQNuROJHat9rYp/17/bxfDmNiUO5HYldZNRZ2Tn72fcNSYAtaz0kD6nfRTQPeuub/ZuPFo8JqRWznZjDuR2JHG3mQ7ioOnPPWTrbgTif1o7G8e7hGt/9f2JKsuHHStdBHJw0T4QOtjjC59O5xOebIjaw1jcicS7zU/R7PxnefXTnmyH+sMpN/NA+1sdLvWTGTtO/X/6pQn+7HSRSTjQPlAdZbz9qcbNYmAgLIf93DW0+maEbsmoOxH/dB84OFLsCcrPUxkL2yEu9btp2tG7JqAshfN/c/WiFCvHHskoOzE4Xn3kV4SiS0C+uu/OAdKy/PakX4Sj5UCenYRia5HLTWTOK0T0KOr8HQ13hRbLwwEWCWgzQ/0OGx4X6ftdEcOz5FKXhbitEpAj8UDle8fbnzcsp8CuiONk55bLwwEWelxdrc74c/3Q/fThkfwAroj3fs2ITYrPlD559ff/rps/GgmG+p+1K8deV2I00oB/fH8ozyU34YNdT967tuEyKwY0OrB9D+/bvdoOxvqftRPfXpdiNOKAb0c77ueP78axkTz0pHXhTit+ZlI5Sd7nI0DzV770rvXhTit9JlIRTLPh+Iq0pafkGRD3YXq7qPn37dcGgi21kD6Wzqvf1zTeT64iJS32m3vj29suTwQbL1bOb897+g0jCljzWfWuY+TqK0T0H/+LKpZPpp+ux1QAd1a+9qRfhK1tR5ndy5Of1a7opuxrW5LLUmLByqzGvkkNQLKShyskx4BZRXySYrWDGhxDWnDK0gXAd2MfJKkFQJ6LAfOVx/rseHjQAV0G/JJohYP6Ln6CI/nY+k9UDkvjt5J1tIBPVefgXQbCnr98/bHdndyCuj65JOELRzQ22nP+6NDyts5qzH1G7Epr00+SdnCAX3e+H6qvrqW1L3wuajl011HJGjhgB7vu53FjufzK09jysOhp59eAVKybECfsbzud1ZPAT16HmgWGrX05CWStGxAf32vDtdrD7ET0Cy0nxry/HKTxYElrBXQ2sfBzxHQ+3Odxp8KsPWupHWs7unzpGnxgN4j9zwF+vzeWLUfrB4sOvpylI13FZ1znQJKmhY/B/oYxFTFL/gq/DOgx+cViZEttvGuoOdSkYCSpoWvwp8O1QdyPh4Degx9IugjoOfqoaK3A/lxMbbxLq/vcruAkqaFA3q7f/NHkbzqCP50CD0FWgX0+udjDqfHfD9cOBvvwg419e/2fw1xW/pWzsfR9n2vs7j4E/hI+iqg59oc/vlz3C6ojXdZ9XA2r8Jfer+GuC3+MJFj7Vzl/ZpeFn8AABHzSURBVHkiofchVQFt7HWO/IxkG++SWjuehjGRvOUfZ1dU88fzy+D7OKuAHusBPY87IWDrXVARz4F9TQPpSdKqT6Sfdhv8M6C1ZgroXpT7nsNnPt3KSXoi+kiP3kN4Ad2HRxoHz3bqJwmKMKDn+vPwjs6B7kDj0lH929ssDqwlwoDWr7yfR17Tt0kvYeDKkX9tkhdXQO/joR7VLJ9yP7woPZZaunwZ8Um+4gtoochmMULq1dPtBXR53XuODFgiHxEF9FJraBHQ05t+9rBNz6vn/5MMWCIjcQX0rnoyyWn8oCjb9Kx6d+nt7ZOPGAM6gY16RkON1E+yIaCE0UiIN6Bjn8N0Z5OfiXzCRUAJIp9wI6CMJp9wJ6CM5OgdKgLKKPIJTwLKGPIJNQLK57q3bdohJWsCyqd6bntXUPImoHyoXUo3vUO0AQ1jUw/V2dH02CUQUD7R+9Sl/q8hIwLKW73nOQUUBJS3+i8TCSgIKG8MP7Su/2vIiIDyyvAoJQEFAeWFV4M8BRQElGEvx8gbxgQCypB3txgZSA8CSq8PbtF0KyfZE1B6NLI4mEn9JHcCSlejinY0YYiA0tZspVOdMEhAaToM9NM/HnQIKHWdQ3XDPWGYgFLjoUswhoDy4KFLMI6AUhoaqtT/NSCglIbGKQkoDBNQbjx0CQIIKC9vezeMCYYJKK/vMjKQHgYJaO7e3qTpVk4YIqCZ+yCM+gkDBDRrughTCGjG7FfCNAKaLfmEqQQ0V/IJkwlonuQTZiCgOXL0DrMQ0PzIJ8xEQLMjnzAXAc1Mfz6NlYcQApqVgUq6WxOCCGhGhp/56XkhEEJA8zG0h+mJdRBIQHMxfIDumckQSEDz8OoEp4BCIAHNwbtHJvd/DbwhoBl4c3ldQCGQgCbv7egkAYVAApq4D0Z3CigEEtCkfTQ43jAmCCSgKfvw3iID6SGMgKbr81sz3coJQQQ0VaN6qJ8QQkDTJIewAgFNknzCGgQ0QfIJ6xDQ5Dh6h7UIaGLkE9YjoGmRT1iRgKZEPmFVApoOR++wMgFNhXzC6gQ0EfIJ6xPQJMgnbEFAE+DoHbYhoNGTT9iKgMZOPmEzAho3+YQNCWjMHL3DpgQ0XvIJGxPQaMknbE1AIyWfsD0BjZKjd9gDAY2QfMI+CGh85BN2QkBjI5+wGwIaF0fvsCMCGhP5hF0R0IjIJ+yLgEZDPmFvBDQSjt5hf+IL6PFQ+TL+h2NtkHzCHkUW0NOhaWxDI62QfMIuRRXQn18Pbb//a9QcouyQfMJOxRTQf/5sFvPX9+vf//h7zCwiLJGjd9itmAJ67uTyltQfY2YRXYrkE3YspoAeuwfs153QUadBY4vRPPmszndMnxNQF1FAr7ub3zrfPI07ho8rIjNF73nGeIaZAU8RBfS6t9k9XD+Pu4wUU0PmSt5jLgoKMxPQnZptj7E2m3jWHuIQUUCzOoSf74C7Pp9oVh/iEFFAM7qINOf5SgGFxcQU0P5hTN290heiKMi813sEFBYTU0B7B9KPuxUpgoLMfblcQGExMQX0Xsym3/4aNYf9F2T20UYCCouJKqD1RzGl+TCRBQZrCigsJrKAXpJ+nN0ig90NY4LFxBfQSXadkIXuFTKQHpYioHux3K2WbuWEhQjoPizaN/2EZQjoHsgbRElAd0A+IU5JB7QzanSXodrnUgHvCejGdrlQwEciCmjPfUg3UT/OTj4hZgK6pUn53Os+NeQjooD2fapx1AGdVr/9npWAbMQU0GIfdNTT6zp2VJuyfaEVfPyIgsJmogroraAjH7/UspvYVNEM3Y+s/cBu1gmyE1dAb0fxoz7Co20vsamKGbwf6RlLsAORBfRyOvR8stzn9tGaxw5n+H6kgMIOxBbQ60H8lF3QPbSmdsAenkEBhR2ILaATd0G3b03jfKeAQtSiC+g0m7emeblIQCFqArrur28ugIBC1AR0zV/evtY+JaCGMcHmBHS9X90dqjRhP9JAethetAE9BQ2p3y42vWPlp+xHupUTNiegqxjq3JT9SP2ErQnoCl50zn4kRExAF/c6j/oJ8RLQpakjJEtAlyWfkDABXZKDc0iagC5HPiFxAroY+YTURRvQMOslTT4hfQK60O/RT0ifgC7yW+QTciCgS/wS+YQsCOj8v0I+IRMCOvcv0E/IhoDOO3v5hIwI6Kxzl0/IiYDOOG/5hLwI6Gxz1k/IjYDONF/5hPwIaO9UIx9yLJ+QJQHtm2hcQeUTMiWgPdOM+qA3R++QLQHtTjLmo4blEzImoK8meTu5fELOBPTVJG8ml0/Im4C+muTl5I7eIXcC+mqSd5/nPsciAdES0FeTDE8un4CAvpxkaHL5BC4C2jfJu2FMjt6BgoD2TPNyIL18AiUB7Zvoxa2c8glUBLR3qqF+yifwJKCjflo/gScBHfGz8gnUCejnPyqfQIOAfvqD8gm0COhnP6afQIeAfvJD8gn0ENAPfkY+gT4C+vYn5BPoJ6BvptdPYIiAvpxaPoFhAvpqYvkEXhDQ4UnlE3hJQIcm1E/gDQHtn0w+gbcEtHcq+QTeyy6gAPOZvVFzz3BOW/9jA2mZvVFzzzAOmRz0Z7KauaxnJqsZ1XpGtKhziuklmiCT1cxlPTNZzajWM6JFnVNML9EEmaxmLuuZyWpGtZ4RLeqcYnqJJshkNXNZz0xWM6r1jGhR5xTTSzRBJquZy3pmsppRrWdEizqnmF6iCTJZzVzWM5PVjGo9I1rUOcX0Ek2QyWrmsp6ZrGZU6xnRos4pppdogkxWM5f1zGQ1o1rPiBZ1TjG9RBNkspq5rGcmqxnVeka0qHOK6SWaIJPVzGU9M1nNqNYzokWdU0wv0QSZrGYu65nJaka1nhEt6pxieokmyGQ1c1nPTFYzqvWMaFHnFNNLNEEmq5nLemaymlGtZ0SLOqeYXqIJMlnNXNYzk9WMaj0jWtQ5xfQSTZDJauaynpmsZlTrGdGizimml2iCTFYzl/XMZDWjWs+IFnVOMb1EE2SymrmsZyarGdV6RrSoc4rpJZogk9XMZT0zWc2o1jOiRZ1TTC/RBJmsZi7rmclqRrWeES0qwL4IKEAgAQUIJKAAgQQUIJCAAgQSUIBAAgoQSEABAgkoQCABBQgkoACBBBQgkIACBBJQgEACChBIQAECCShAIAEFCCSgAIEEFCCQgAIEElCAQAIKEEhAAQIJKEAgAQUIlG1Af349XP3219bLsYJf33//19bLsKBf32+v5LetF2MNib+Shcg2zEwDet/obtJ/R/7zZ9IreSxfyGg2uXCJv5I30W2YeQb0+TJF80KFOya9jqd8XsjEX8mb+DbMPAN6qo75brsvX7ZemkVd91pieS+GuB3w/fH35XJO/oVM/ZUsxLdhZhnQ2//P/bh/eUr8PXk/o5TuKh6r7ey6omkfxKf+St5EuGFmGdDz8//dbv+3/mPThVlScUT07xM+c3ZdwWrdTklfSEr+lSxEuGFmGdD6tnaM43UKczsi+pHypYfrbtmX55e3Y/lEJf9KFiLcMLMMaN0x5SO/0y0vKW92tS0u5dXM4JXsiGTDzD2g10OjhPdbCilvdqfnbsp1NaPY4iZI+ZVsi2XDzD2gaZ85K6S82dWP8yLZZZkg5VeyLZYNM/OApn3i7C7lza4eTQFNSDQbZt4Bvb5M6b8jU97smgGN4qrDBCm/kk3xbJhZB/QcyVizaVLe7OyBJimiDTObgD5uEntuZRG9TB/rWc2UNzsBTVFMG2bGAT0eIjnNMkZmAXUVPkFRbZjZBvR2p0MUN9uOk1lAz7mMAy1ksIrRbZjZBLTldmdxFMMkpkt5s8vmTqRCyq9kJbYNM9OAJv/oiZqUN7ts7oUvpPxKlqLbMPMM6HW7i+tlmiLpzS6fpzFdEn8lC/FtmFkGNJYnvcwj6c0uo+eBJv5K3kS4YWYZ0POhIbKXbKy0N7ucnkif9it5iXLDzDGgxaO9I3udJkh8s/OZSMmIccPMMaD1D16J5XWaIPXNLp9P5czjlYxrw8wxoACzEFCAQAIKEEhAAQIJKEAgAQUIJKAAgQQUIJCAAgQSUIBAAgoQSEABAgkoQCABBQgkoACBBBQgkIACBBJQgEACChBIQAECCShAIAEFCCSgAIEEFCCQgAIEElCAQAIKEEhAAQIJKEAgAQUIJKAAgQQUIJCAAgQSUIBAAgoQSEABAgkoQCABBQgkoACBBBQgkIACBBJQgEACChBIQAECCShAIAEFCCSgAIEEFCCQgAIEElCAQAIKEEhAAQIJKEAgAQUIJKCMczq0fbv88+fht7/m/TW//sftf68z/v1fH0w9MN3Pr8UCfjSLEYvz6VKRPgFlnHUCerzO9TIxoPUl/THn4ggoFQFlnFUCej5MDuj1W3VfZlwcAaUioAT59b0WkcUC+qlu0u79vO94HqcWtLU4AkpFQAmy94AeD7VZFDUdNcOXiyOgVASUIDsP6O3yUf2857WnE5ZQQBkgoATpCej5dqT8x9+PSc7tyzfNbxyv096+cw9b/b+dH1d+aqlqz704Lq+q2EnasXXMfl3cooG11F8bW82uZ16n56/rLE57qWrreJ5+wpWYCChBugE9Ni94l0OInoOIbg2rf+Ma0P9Tdas5cU9AHxOUc39eyioi1w7o9Xe1djhP9yn7Atozr/9Z/3XDAW2v4+PvU84XEBMBJUgnoP/h0IjJsyXlNx79rL5xPPy7r+XuWmvibkCfP3zPX30oQLlj2QjoubGzWi7P7Ud7Atozr6fbxIMBba9jbRWnjpsiEgJKkHZAy2acahUqSnU+lEe01VWdU/0b98p1J24PYzreZ38r1JfL84C8+KLcsWwE9NQ5jL79jh99Ae2bV2tlBoYxdRb7VIb09h9a/SZRAkqQTkDLv5Xpeu4ClkfT5R7g5fnV8TFJZ+J2sR5H5OUX50cfrzPrG5l56hxEDwa0b16tlRkIaGexj9WOp6tM2RBQgnQCWgarzNLxeQ7yXp/aPmFZt+PjZzoTt4v1jNypdXBcLka7WMfuWcjjQED75tVamYGA9q2jHc/MCChBhoYx3Ztz/a+PlFy/86W2d1Z94/md7sTtYrWzWfu9h6l7oH3zak/SG9DuYp+d/syOgBLkbUDr/vi7PlK0vZPambhTrGP7mvrzys5gQLvnQAcuIvXM69OAthe7GojgInw2BJQgrwNauz79COhj8nLX7Vg7KTo2oPUf6Qto71X4gTr2zOujgPYsdu2KvoTmQUAJ8jagzYC92gPtTPw2oNXO3+//GjgH2jsO9EvfkvbO69OA9pzxfOyXKmgWBJQgbw/hu/emV+cGy0tC9UP49jXrN+dAT4/r5AMBrd2J9Ot7NfDpR9+S9s7r00P4gUvtRURdTsqCgBLkdUDLSzY1tbOSx+oq/PPeyfaFl8Gr8MXVmtqvO/cfwhfH19+qKa4/8hgz9dyZvV8z753XRwHtLHZtIYbjSloElCCvA1rbryvr0zcOtJay5sTtYj1/+Nh8/OhtX6//CZ3PpzGV5yWfM/jy+NFGQGvz+iig3cV+Dp4S0FwIKEHeBLQ4ii32z6rnINXvRCqq9Axod+Jz6x738k6k25/3c6L3b5+qMg49D/Tes+Ozn4+7ix7PJumbV09A64vzHN7fXOzHjVTVHVMkT0AJ8iagtaFBtTsla9e62+PQmxPfr3D33Qv/rF+lP6DNW9oP9ZuGSv/2vbrg35lX/4X6zr3wncU+Pv9uBzQPAkqQdwHtPpmo+zSm54Xy9sTV/uMjjK2L249S/efH+PhusRppPNTPiRZ/qwbC98yrvTKtxek+jamx2FWJyYCAEuRtQKuTj7VD2fbzQOuVaU1cpOhL58mbj195T9e36vrSwM3nZeCuv6d+Yfz4+FY99fV5dVamuTi139ZexzL0bkbKhYCSh2sMZY25CShAIAEFCCSgAIEEFCCQgAIEElCAQAIKEEhAAQIJKEAgAQUIJKAAgQQUIJCAAgQSUIBAAgoQSEABAgkoQCABBQgkoACBBBQgkIACBBJQgEACChBIQAECCShAIAEFCCSgAIEEFCCQgAIEElCAQAIKEEhAAQIJKEAgAQUIJKAAgQQUIJCAAgQSUIBAAgoQSEABAgkoQCABBQgkoACBBBQgkIACBBJQgEACChBIQAEC/X8Gt2eYB66HMgAAAABJRU5ErkJggg==)
# fit linear regression model, see the difference and similarities
fit2<-lm(nitrogen ~ influent.f)
anova(fit2) #compare to summary(fit1)
## Analysis of Variance Table
##
## Response: nitrogen
## Df Sum Sq Mean Sq F value Pr(>F)
## influent.f 5 1925.2 385.04 9.0441 2.202e-05 ***
## Residuals 31 1319.8 42.57
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# fit a random effects model, assuming two nested errors, defaut is REML
fit3 <- lme(nitrogen~1,random=~1|influent.f)
summary(fit3)
## Linear mixed-effects model fit by REML
## Data: NULL
## AIC BIC logLik
## 258.3511 263.1017 -126.1756
##
## Random effects:
## Formula: ~1 | influent.f
## (Intercept) Residual
## StdDev: 7.957598 6.531317
##
## Fixed effects: nitrogen ~ 1
## Value Std.Error DF t-value p-value
## (Intercept) 21.22312 3.429034 31 6.189241 0
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -1.91443594 -0.53645956 -0.03217348 0.83713515 1.95823685
##
## Number of Observations: 37
## Number of Groups: 6
intervals(fit3)
## Approximate 95% confidence intervals
##
## Fixed effects:
## lower est. upper
## (Intercept) 14.22956 21.22312 28.21668
## attr(,"label")
## [1] "Fixed effects:"
##
## Random Effects:
## Level: influent.f
## lower est. upper
## sd((Intercept)) 3.951429 7.957598 16.02544
##
## Within-group standard error:
## lower est. upper
## 5.089555 6.531317 8.381499
fixef(fit3) #overall mean, intercept
## (Intercept)
## 21.22312
ranef(fit3) #alpha_i hat
## (Intercept)
## 1 0.309286
## 2 -6.719332
## 3 -3.897945
## 4 2.946104
## 5 -6.012984
## 6 13.374871
yhat<-predict(fit3)
resid<-resid(fit3)
plot(resid~yhat)
![](data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAABUAAAAPACAMAAADDuCPrAAAATlBMVEUAAAAAADoAAGYAOpAAZpAAZrY6AAA6OgA6ZmY6kLY6kNtmAABmtv+QOgCQkNuQ29uQ2/+2ZgC2///bkDrb////tmb/25D//7b//9v///+Xg4BjAAAACXBIWXMAAB2HAAAdhwGP5fFlAAAfRElEQVR4nO3d62KbVsCmUbWd9JtJ05mqrpv6/m90dLSRRYh4QbA3e60/cVxChRGPOWzQ7g2AyG7tFwBQKwEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQKjsgO4AZjN/omaf44zW/mkD2zJ7o+ae4Zye8AsDaJaAAoQEFCAkoAAhAQUICShASEABQgIKEBJQgJCAAoQEFCAkoAAhAQUICShASEDjOT3nYVZAPQQ0nZGCQvMENJzP7tMXQHsENJvN7uPLWeYIVEhAp85GQaFZAjp1NgIKzRLQqbMRUGiWgE6djYBCswR06mwEFJoloFNnI6DQLAHNZmMYEyCg6XwMpAcENJ2RWzmheQIaz0k/oXUCChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCNUX0P3u6sv4fyygwHwqC+jL7tbYhgooMJ+qAvrv77vPfv171BwEFJhPTQH978/bYn7/4/D33/4ZMwsBBeZTU0Bf73J5TOq3MbMQUGA+NQV0f3/AftgJHXUaVECB+VQU0MPu5te7b76MO4YXUGA+FQX0sLd5f7j+Ou4ykoAC8xFQgFBFAXUID5SlooC6iASUpaaA9g9jut8rHSCgwHxqCmjvQPpxtyIJKDCfmgJ6LuatX/4aNQcBBeZTVUC7j2LyMBFgbZUF9M3j7IBi1BfQEe4O+HcCCsxHQAFCmw7oPQEF5iOgACEBBQgJKEBIQAFCFQW05z6k0R8rJ6DAfAQUIFRRQPs+1VhAgfXUFNDTPuiop9fdEVBgPlUF9FjQkY9f+kRAgfnUFdDjUfyoj/D4TECB+VQW0LeXXc8nyz1OQIH51BbQw0H8lF1QAQXmU1tAJ+6CCigwn+oCOo2AAvMRUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKC9U509+9UAdRPQvokUFHiAgPZMs/v0BUAfAb2fZPfx5VNfDFA5AR2aREGBAQI6NImAAgMEdGgSAQUGCOjQJAIKDBDQoUkEFBggoEOTCCgwQEDvJzGMCXiIgPZMYyA98AgB7ZvIrZzAAwS0dyr9BH5OQAFCAgoQElCAkIAChAQUICSgACEBBQgJKEBIQAFCAspK3O9F/QSUdXjiAIUa88YUUFbhmVcUatSvdgFlDZ66SqHG/WoXUNbguf+UaeSvdgFlDQJKmUa+MwWUNQgoZRLQIbbVQggoZRLQIbbVQggoZRLQIbbVQggoZRLQIbbVQhjGRJkEdIiNtRQG0lMkw5iG2FiL4VZOimQg/QBbazn0kyK5lfPHbK7cknE+8zCRH7KhcMOJBCYRUBrmUhbTCCjtMpiKiQSUdhnOz0QCSrsElIkElHYJKBMJaO9ULsw2QUCZSED7JlLQNggoEwlozzSGtjRCQOlhIP0PPfR4FUNbWmFdc8+tnD828vlUtqqNc7TBZx4mMkBAueF8N7c8zm6IgHJLP7nhgcpDBBQYIKBDBBQYIKBDBBQYIKBDDGMCBgjoEAPpgQECOsStnMAAw5iGeJgIMMRA+gGiCAxyK+ePCSi3HG3wmYeJ/JANhRvOdzOJgNIwIy6YRkBplzG/TCSgtMtdZ0wkoLRLQJlIQGmXgDKRgNIuAWUiAWUlBYwfElAmElDWUcIITAFlIgFlFUWMwDSMiYkElDUUkq4iMk7FBJQ1lHLwXMKJBComoKyhlICWcCmLigkoaygmoDCFgLIGAWUTBJQ1CCibUGlAv/9xPG/12z9j/51ttRACyiZUFNBDNK/F3F9P/X8ZOQ/baiEKGcYE01QZ0P3H4JORO6E21lIYgckW1BjQ18M29/XynZH7oDbWYhiByQZUGNDDn7/+ffney+6Xv8bMw9ZaDv2kfhUG9PWy/3n035/jdkFtrsB8KgzozV7ny7izoAIKzKfCgO67AX39OJzveSk9nvXqgPZUGdBOMwUUWE2FAX15eA/0noAC86kwoK+73bf3b+6dAwXWUmFAu1feu1fkHyGgwHzqCujuNIT+vZqHlI47ghdQYEb1BfTklM3THZ3ffvrvugQUmE9FAX3rNPQU0JfR/RRQYEZ1BfTs39+vl+PHPoxJQIEZ1RjQCQQUmI+AAoQEFCAkoAAhAQUICShASEABQgIKEBJQgJCAAoQEFCAkoAAhAQUICShASEABQgIKEBJQgJCAAoQEFCAkoL1TnT371QB1E9C+iRQUeICA9kyz+/QFQB8BvZ9k9/HlU18MUDkBHZpEQYEBAjo0iYACAwR0aBIBBQYI6NAkAgoMENChSQQUGCCgQ5MIKDBAQO8nMYwJeIiA9kxjID3wCAHtm8itnMADBLR3Kv0Efk5AAUICStMcbTCFgMZzsuFtgPPdTCKg6YxseBtgxAXTCGg4nyU3PK1+EmN+mUhAs9ksueHZ230Wd50xkYBOnc3TNzyHmU8joEwkoFNn8+wNz2Hm8wgoEwno1Nk8P6DL/b9a42fLRAI6dTYCWi8/WyYS0KmzEdB6+dkykYBOnY2A1sv5ZSYS0Gw2y214AvpERjgwjYCG81lswxPQZzLGlkkENJ3RUhuew8yn0k+mENB4TktteA4zoVQCWj6HmVAoAa2AfkKZBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBLQCBtJDmQS0fG7lhEIJaPE8TARKJaCl8zg7KJaAls4DlaFYAlo6AYViCWjpBBSKJaClE1AoloCWTkChWM8N6Pc/dn1+/Xvu/+nDL66+AgkoFEtAS2cYExRLQItnID2UaplzoK+73S9/nb46JvXb3P/Lx1WZILdyQqEWCei/v+++vv/l9drSNdTZIP2EMi0S0P3uyw//tiwRAuazREAPh+3do/ZX50CBTRBQgNASAf3vz0+H8L/9M/f/9FECCsxnoXOgnV3Ql51zoMAmLHUV/nrYfvzSVXhgE5YZB/pyM47eOFBgExZ6mMjrRz5X3P8UUGBOiz2N6XJX53oX4E8EFJiPx9kBhAQUICSgAKHnP87u2/1D7dyJBGyBgAKEBBQg5BwoQMeY5+8KKMCHUZ8AIaAA78Z9BtmyAT0+SmS9Z9m9CSgwaOSn4C4U0O9/HLt57Oe6N8MLKDCgm4gHcrHYp3L++vfxwcpr3w8voMCAEgN6eQjo4Y9DO192nY/oXJqAAgNKDOjlIfSXdE78VM79+zWyYDYCCgwoMKCHQ/fjGdDDH6fTn1M+VO72yczjGyqgwIACA3r5VM7DH6dy5gE9X4SaclOTgAIDyg3o62WHMQ7o6SJU59+ebhEdNypKQIEBBQ5jugT0+tmcL+lQ0Ne7XB6TOuoTlgQUGFLeQPpD5r5+nAI95DS8iLS/33UdOzMB5daY+55pQnm3cu6P6Xw57z8ez2NmH8t57vAnI3dn51tcG94mjNpYaENxDxO5Xvz5dh6EFB7BX84E3Bo8oXp3yWm+7cSGtwnjDtfgsyU/F/4Yzpf8RqSSAmrD24SRFwzgswU/F/5yCT47fn8r6hDehrcNI4eswGc1Pc6unItINrxtsB6ZqKaA9g9jGnVjvYDSYT0y0WIBPY56P16Ln/A40N6B9OPOqAooHdYjEy32PNDjBZzzYKb4JOj9p9ONfriogJajgGEM1iMTLRPQS/kOtdvvphS08yims5UeJmLDm66EgWDWIxMtEtDjsfdv//z7+3F3cT/1gcolPM7OhjdZEQPBjKagR3ED6V9Pl3rOAd3EA5VteFMV8hMsIuOUZdSx0UK3cp5v4jwF9PJ00HUYSF+IUvbhSziRQFHGbdvLPUzkGtDDLuh6H4rkVs5ClBLQEi5lUZKRx0YLPs7uGtApT6SfysNEClFMQOFGuQ9U3lZAmURAKVOBAf10CL/fwDlQJhJQylRgQC83sV8Cmj9QeQa21UIIKGUqMaDnT0M6B3T0p3DMyrZaiEKGMcEnJQb0PHr+FNDjk0HXO4K3sRbDQDCKVGRAb25iX+8SkoAWxEAwSlTgMKa3y4OU1u+ngBZEPylReQPpTy4JXTWfAgr8RHm3cpZDQIFhxT1MZL/yfucHAQXms9BA+hUvvN8QUGA+C93Kud7Q+VsCCszHHihAaKk7kda7+eiGgALzWeYq/L+/F7IPKqDAfBY6B3qr9MfZGeINPEJA+yZS0AX4IVM/Ae2ZxmMuFuDXFBvgTqT7STxobQF+TbEFAjo0iW37WfyaYhMEdGgS2/az+CmzCQI6NIlN+1n8lNkEAR2axKb9LH7KbIKADk1i034WP2U2QUCHJrFpP4ufMpsgoPeTuD78fALKJghozzRGKD6dX1NsgoD2TeQemafza4otENDeqfTz6fyaYgMElJXoJ/UTUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUFZiID31E1DW4VZONkBAWYWHibAFAsoaPM6OTRBQ1uCBymyCgLIGAWUTBJQ1CCibIKCsQUDZBAFlDQLKJggoaxBQNkFAWYNhTGyCgLIKA+nZAgFlHW7lZAMElJXoJ/UTUJom40whoLTMiQQmEVAa5lIW0whoBewkPYnBVEwkoOVzmPkshvMzkYAWz2Hm0wgoEwlo6RxmPo+AMpGA9k5V0BGzjfx5/GyZSED7JiqpoDby5/GzZSIB7ZmmqHOONvLn8bNlIgG9n6Ssc4428ucpbV1THQEdmqSEraq017MphR1tUB0BHZqkhK2qtNezLWWd76Y6Ajo0SQmblcPMp9JPphDQoUmK2K4cZkKpBHRokjKK5TATCiWgQ5MUkiz9hDIJ6P0kzjkCDxHQnmmccwQeIaB9EzlmBh4goL1T6SfwcwIKEBJQgJCAAoQEFCAkoAAhAQUICShASEABQgIKEBJQVuJ+L+onoKzDEwfYAAFlFZ55RanG/GYXUNbgqauUatSxkYCyhgKf+w9H446NKg3o9z+OvyB++2fsv7OtFkJAKdPIY6OKAnqI5rWY++s+9peR87CtFkJAKdPId2aVAd1/nKUYuRNqWy2EgFKmBgL6egjn18t3Ru6D2lYLIaCUafsBPfz569+X773sfvlrzDxsq4UQUMq0/YC+XvY/j/77c9wuqG21EIYxUabtB/Rmr/Nl3FlQG2spDKSnSNsP6L4b0NePw/lH2FiL4VZOSrT9YUz7bjMFtFr6SYk2O5C+9xBeQIE5bfVWzs5FpG/v39w7BwrMaaMPE7kGtHvlvXtF/hECCsynroDuTkPo36t5SOngEfyux7NeHdCe+gJ6csrm6Y7ObwP/QkCBZ6oooG+dhp4C+vKTfvYQUGA+dQX07N/fr5fjxz6MSUCBGdUY0AkEFJiPgPZO5Xwp8HPVBnTsc5jOHhwbq6DAAwS0ZxqPuQAeIaD3k3jQGvAQAR2aREGBAQI6NImAAgMEdGgSAQUGCOjQJAIKDBDQoUkEFBggoEOTCCgwoNqAZgxjAuYjoD3TGEgPPEJA+yZyKyfwAAHtnUo/gZ8TUICQgAKEBBQgJKAAIQGtgItaUCYBLZ9hVVAoAS2egf1QKgEtnVtLoVgCWjoPN4FiCWjpBBSKJaClE1AoloCWTkChWAJaOgGFYglo6QQUiiWgpTOMCYoloMUzkB5KJaDlcysnFEpAK6CfUCYBBQgJKEBIQAFCAgoQElCAkIAChAQUICSgACEBBQgJKEBIQAFCAgoQElCAkIAChAQUICSgACEBBQgJKEBIQAFCAgoQElCAkIAChAQUICSgACEBBQgJKEBIQAFCAgoQEtB4TmdzzQ+oj4CmM1JQaJ6AhvPZffoCaI+AZrPZfXw5yxyBCgno1NkoKDRLQKfORkChWQI6dTYCCs0S0KmzEVBoloBOnY2AQrMEdOpsBBSaJaDZbAxjAgQ0nY+B9ICApjNyKyc0T0DjOekntE5AAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUIBQfQHdvz/G48v4fyygwHwqC+jL7tbYhgooMJ+qAvrv77vPfv171BwEFJhPTQH978/bYn7/4/D33/4ZMwsBBeZTU0Bf73J5TOq3MbMQ0HJ4oir1qymg+/sD9sNO6KjToDbXYnimPxtQUUAPu5tf7775Mu4Y3tZaCp8qxRZUFNDD3ub94frruMtINtZC+FxTNkFAWUN3RVgpVKuigDqE3xABZRMqCqiLSBsioGxCTQHtH8Z0v1c6wLZaCAFlE2oKaO9A+nG3ItlWCyGgbEJNAT0X89Yvf42ag221EALKJlQV0O6jmDxMpGqGMbEJlQX0zePsNsJAeragvoBOYmMthls52QABZSX6Sf0ElKbJOFMIKC1zIoFJBJSGuZTFNJsO6N2oUZsJXQZTMZGA0i7D+ZmoooD23Ie0G3svp62EDgFlIgGlXQLKRBUFtO9TjQWUCQSUiWoK6GkfdNTT6+7YSugQUCaqKqDHgo58/NInthI6BJSJ6gro8Sh+1Ed4fGYrocMwJiaqLKBvL7ueT5Z7nM2ELgPpmaa2gB4O4qfsgtpMuGGAMJPUFtCJu6C2E27pJ1NUF9BpbCjAfAQUICSgACEBBQhVG9CXaEi9gALzEVCAkIAChAQUICSgACEBDbmDBRDQjHuoAQHNeIoPIKAZz5EE3ioOaGa2gM4+S6A+Ajp1NgIKzRLQqbMRUGiWgE6djYBCswR06mwEFJoloFNnI6DQLAHNZmMYEyCg6XwMpAcENJ2RWzmheQIaz0k/oXUCChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAe6fyrE/g5wS0byIFBR4goD3T+Lwj4BECej+JT9wEHiKgQ5MoKDBAQIcmEVBggIAOTSKgwAABHZpEQIEBAjo0iYACAwR0aBIBBQYI6P0khjEBDxHQnmkMpAceIaB9E7mVE3iAgPZOpZ/AzwkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAoeYCCjCf2Rs19wzntPYPG9iW2Rs19wzL0vhJAIvfsrYXf6Gl3/jPuO33kMVf+xWsqu3FF9A5tP0esvhrv4JVtb34AjqHtt9DFn/tV7CqthdfQOfQ9nvI4q/9ClbV9uIL6Bzafg9Z/LVfwaraXnwBnUPb7yGLv/YrWFXbiy+gc2j7PWTx134Fq2p78QV0Dm2/hyz+2q9gVW0vvoDOoe33kMVf+xWsqu3FF9A5tP0esvhrv4JVtb34AjqHtt9DFn/tV7CqthdfQOfQ9nvI4q/9ClbV9uIL6Bzafg9Z/LVfwaraXnwBnUPb7yGLv/YrWFXbiy+gc2j7PWTx134Fq2p78QV0Dm2/hyz+2q9gVW0vvoDOoe33kMVf+xWsqu3FF9A5tP0esvhrv4JVtb34AgpQNgEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQGijAf3+x69/X7/+78/du29rvqhl/Pv7cUF/+evjO9//OH7n63ovaUmfF7/Jtd9d1y2t/c9Lv8DK32ZADz+4j4Ce30GNbEIfC/v+A9hfvtFt6lbdL35Ta/8jGO+L2tDav1/6BVb+NgO67/Tj7XX39J9iMbrvmMtP4OXzNzasZ/FbWvs9O1wNrf2epV9g5W8xoKef5Mf75WX7W8674/ZyOoI57nh8OX5xPKr57Z/ze+nLui/u+e4Xv8m1f/w9clznba39+6VfYuVvMKDnEyEfAd1v/5fv1fG987HrcVrs/XXLOfxYtn4Y17P4La39w47DZQ0ffhDnrxpa+z1Lv8TK31xAT0dx/6tzDvTwg738Ptq+14/9jONu+LfTj+P6k3jZ/KWE+8Vvau0fInld/MuuV0tr/37pF1n5mwvocUf+W/ci0uFNtPWDl3fdrWR/ehd13laHLzfekvvFb2rtd1yWvqm1/+Gy9Ius/A0G9MvtVfjDbsm3027pxg9hPtufFrjTlJuhCZt3XvxG1/51z7PNtf++373Eyt9cQE+675aX3S//c7kQ18yv4LfTu+i4uJ3T6B8niRpwWfw21/7L9aJzk2v/fekXWfnbD+h1IFxb29B152N/MySwkU3o7WPfq721v+8M2mlv7XeXfpGVv/mAHq8mnN85r7utn0fvuJ7y6m42bWxCJ9fFb3Dt3wa0tbXfXfpFVv7mA/oxpuG4VbVyGuh9UW83oVaGRL4vfntr/zqcvOfXZwNr/2bpF1n5mw9oVxPvoaPX93Gw7e2DvHUXv6uZtX8ey3e8/tzk2n9f+q6nrfymArr5sXAXnYC0uAn197OZtX902flqce2/3ex6Xj1t5TcV0Nc2NqF955R5g9dh9z+4YtDI2j87B6PBtX9yn8unrXwB3ZrjaaCPA5jX1kYC3i5+VxNr/+q8sM2t/Yv7VS2go/z4EH7zZ8GOTwL4evPXpu5F+bT4XS2s/XfnYLS29q/uc/m0lb/5gN7ejLH5o5jPD41o6W7ot/vFb2rt39/J2tLav1/6RVb+5gPaGcDwuv1Het2fPm/oeTw9i9/U2u8s7PXLhtb+/dIvsvI3H9DjObHz1y/bvyH6+gyijpaeCHm/+M2t/feFPe98NbT275d+kZW/+YBeHg96tu2DmNsncHfuh77Y/EWEnsVvae3fLOwllw2t/fulX2Llbz+gnR/j1q8hdD/UoHs/39nWd8D6F7+dtf/W/UiT993NZtZ+39IvsPIbCOh1z2TjhzBvtx8J1HnXtPK5jD9Y/FbW/sk5Gd3r7a2s/aP7pX/6yt9mQAEWIKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUKr035+7X/766VTf/+8CL4WGCShVeiig+93XJV4L7RJQqvRIQF93AspzCShVElBKIKBUSUApgYBSJQGlBAJKFTox/P7H4ctzQA/f3e1+++d9qv3x75eynv7bwbc1Xi6NEFCqcKjmNZSvx0KeArq/beTL7uo4qYDyfAJKHfbXQ/ZDOg99PAb0f669/PXv43/46Odpb1VAeT4BpQ7//n45hj9/cQjoJY4v514e91HPUxy/OO2tOgfKswkodTjveL4dg3ncFT0G9LzjeSzol7dTLr+cJz0k9vSfBJRnE1AqcQ7nNaTHgF7qeOhl5zrS22kXVEBZhIBSicsx/OWPzjCmTwE9HdwLKIsQUCpx2fV8OcexN6Cvu5vLSgLKswkotbgOXzqd6ewJ6OFPAWVZAkotTgPoD5U8jUu6D+jx4vu5nc6BshQBpRr7QygvR/A9AX15vywvoCxFQKnG4Rj+/12O4O8D2vnGq0N4FiKgVOOwZ/m/L0fwQwE9HssLKIsQUOrx8n7bZs8h/P7y317enyfyuvs0QBRmJqDU43id/XK30X1AP8YwXQN6vizvXnieR0CpR6eaPcOYrs9m2v2fP8/ZPN8v7yie5xFQKrJ/PybvG0h/3uP8+nFb/KmgX9Z5qTRBQKnH+xOZoAwCSj1eHvkseFiOgFKNww6o43GKIqDU4dXj5SmPgFKHU0DtgFIWAaUOx0vsriBRGAEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUL/H2MLN6bQklQ/AAAAAElFTkSuQmCC)
par(mfrow=c(2,2))
plot(fit3) #standardized residuals v.s fitted values
![](data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAABUAAAAPACAMAAADDuCPrAAAAbFBMVEUAAAAAADoAAGYAOmYAOpAAZmYAZrYAgP86AAA6OpA6kJA6kLY6kNtmAABmAGZmOgBmtrZmtv+QOgCQkNuQtpCQ29uQ2/+2ZgC2/7a2///bkDrb/7bb/9vb///m5ub/tmb/25D//7b//9v///8X83cZAAAACXBIWXMAAB2HAAAdhwGP5fFlAAAgAElEQVR4nO3dbWPb1podUOZO6nhae6Z226hX01zZ0f//j+W7KFncEiHg4ADPWl8sOYkRbPrZAgngYPMIwCCbuf8HAJZKgQIMpEABBlKgAAMpUICBFCjAQAoUYCAFCjCQAgUYSIECDKRAAQZSoAADKVCAgRQowEAKFGAgBQowkAIFGEiBAgykQAEGUqAAAylQgIEUKMBAChRgIAUKMJACBRhIgQIMpEABBlKgAAMpUICBFCjAQAoUYCAFCjCQAgUYSIECDKRAAQZSoAADKVCAgRQowEAKFGAgBQowkAIFGEiBAgykQAEGUqAAAylQgIEUKMBAChRgIAUKMJACBRhIgQIMpEABBlKgAAMpUICBFCjAQAoUYCAFCjCQAgUYSIECDKRAAQZSoAADKVCAgRQowEAKFGAgBQowkAIFGEiBAgykQAEGUqAAAylQgIEUKMBAChRgIAUKMJACBRhIgQIMpEABBlKgAAMpUICBVlygmxXvW2THyym757Pv+Nzbn9Ds2c7FjpdTds9n3/G5tz+h2bOdix0vp+yez77jc29/QrNnOxc7Xk7ZPZ99x+fe/oRmz3Yudrycsns++47Pvf0JzZ7tXOx4OWX3fPYdn3v7E5o927nY8XLK7vnsOz739ic0e7ZzsePllN3z2Xd87u1PaPZs52LHyym757Pv+Nzbn9Ds2c7FjpdTds9n3/G5tz+h2bOdix0vp+yez77jc29/QrNnOxc7Xk7ZPZ99x+fe/oRmz3Yudrycsns++47Pvf0JzZ7tXOx4OWX3fPYdn3f7G4Cmxq2wUf+0G7cN0Ni4JTbqn3bjtuc+/O7Zv/419/9Bz6STSCdQoDUYgkQ6iXQCBVqDIUikk0gnUKA1GIJEOol0AgVagyFIpJNIJ1CgNRiCRDqJdAIFWoMhSKSTSCdQoDUYgkQ6iXQCBVqDIUikk0gnUKA1GIJEOol0AgVagyFIpJNIJ1CgNRiCRDqJdAIFWoMhSKSTSCdQoDUYgkQ6iXQCBVqDIUikk0gnUKA1GIJEOol0AgVagyFIpJNIJ1CgNRiCRDqJdAIFWoMhSKSTSCdQoDUYgkQ6iXQCBVqDIUikk0gnUKA1GIJEOol0AgVagyFIpJNIJ1CgNRiCRDqJdAIFWoMhSKSTSCdQoDUYgkQ6iXQCBVqDIUikk0gnUKA1GIJEOol0AgVagyFIpJNIJ1hCgd5vtj69vW0Fep0hSKSTSCfov0B/ft3s/fbnW9tWoNcZgkQ6iXSC7gv07+/76vzxefP7X29sW4FeZwgS6STSCbov0Ifjoee2Qb+9sW0Fep0hSKSTSCfovkDvjkee2yPRNz4GVaCBIUikk0gnWECBfnrxxdVtK9DrDEEinUQ6QfcFevLz6+bLG9tWoNcZgkQ6iXSCxRTowyun4TcvTLbx5TMEiXQS6QQvO+iDJTRZh20PQH89C69A380QJNJJpBMspEC3/fmPf7617TYF+sdBk22NxxAk0kmkEyzjLfy2P9+8jr5Rgf7xxyIb1BAk0kmkEyyiQHc3I71xEehjowI9NefSGtQQJNJJpBMsoUB/fH7H8WezAn35xTIYgkQ6iXSCBRToO/uzSYFe1OayGtQQJNJJpBP0X6Db/nzz/NFh2wr0OkOQSCeRTtB9gb7n/Ptx29MX6GVpKtD1kE4inaD7Ar1/3/v3R0egkSFIpJNIJ+i9QE+rgb5jTWUFGhiCRDqJdILeC/TH564K1Fn4VZJOIp2g9wK9ZduuA73OECTSSaQTKNAbuRNphaSTSCdQoLdaZH8agkg6iXQCBVqDIUikk0gnUKA1GIJEOol0AgVagyFIpJNIJ1CgNRiCRDqJdAIFWoMhSKSTSCdQoDUYgkQ6iXQCBVqDIUikk0gnUKA1GIJEOol0AgVagyFIpJNIJ1CgNRiCRDqJdAIFWoMhSKSTSCdQoDUYgkQ6iXQCBVqDIUikk0gnUKA1GIJEOol0AgVagyFIpJNIJ1CgNRiCRDqJdAIFWoMhSKSTSCdQoDUYgkQ6iXQCBVqDIUikk0gnUKA1GIJEOol0AgVagyFIpJNIJ1CgNRiCRDqJdAIFWoMhSKSTSCdQoDUYgkQ6iXQCBVqDIUikk0gnUKA1GIJEOol0AgVagyFIpJNIJ1CgNRiCRDqJdAIFWoMhSKSTSCdQoDUYgkQ6iXQCBVqDIUikk0gnUKA1GIJEOol0AgVagyFIpJNIJ1CgNRiCRDqJdAIFWoMhSKSTSCdQoDUYgkQ6iXQCBVqDIUikk0gnUKA1GIJEOol0AgVagyFIpJNIJ1CgNRiCRDqJdAIFWoMhSKSTSCdQoDUYgkQ6iXQCBVqDIUikk0gnUKA1GIJEOol0AgVagyFIpJNIJ1CgNRiCRDqJdAIFWoMhSKSTSCdQoDUYgkQ6iXQCBVqDIUikk0gnUKA1GIJEOol0AgVagyFIpJNIJ1CgNRiCRDqJdAIFWoMhSKSTSCdQoDUYgkQ6iXQCBVqDIUikk0gnUKA1GIJEOol0AgVagyFIpJNIJ1CgNRiCRDqJdAIFWoMhSKSTSCdQoDUYgkQ6iXQCBVqDIUikk0gnUKA1GIJEOol0AgVagyFIpJNIJ1CgNRiCRDqJdAIFWoMhSKSTSCdQoDUYgkQ6iXQCBVqDIUikk0gnUKA1GIJEOol0AgVagyFIpJNIJ1CgNRiCRDqJdAIFWoMhSKSTSCdQoDUYgkQ6iXQCBVqDIUikk0gnUKA1GIJEOol0AgVagyFIpJNIJ1CgNRiCRDqJdAIFWoMhSKSTSCdQoDUYgkQ6iXQCBVqDIUikk0gnUKA1GIJEOol0AgVagyFIpJNIJ1CgNRiCRDqJdAIFWoMhSKSTSCdQoDUYgkQ6iXQCBVqDIUikk0gnUKA1GIJEOol0AgVagyFIpJNIJ1CgNRiCRDqJdAIFWoMhSKSTSCdQoDUYgkQ6iXQCBVqDIUikk0gnUKA1GIJEOol0AgVagyFIpJNIJ1CgNRiCRDqJdAIFWoMhSKSTSCdQoM/8cfDxP6g3hiCRTiKdQIFe+uOPtTaoIUikk0gnUKAXTs25wgY1BIl0EukECvTCuTcVaC3SSaQTKNAnF7W5ugY1BIl0EukECvSJAq1KOol0AgV6dlmaCrQU6STSCRToE0egVUknkU6gQJ8o0Kqkk0gnUKAXnIUvSjqJdAIFesF1oEVJJ5FOoEAvuROpJukk0gkU6DNr7U9DEEknkU6gQGswBIl0EukECrQGQ5BIJ5FOoEBrMASJdBLpBAq0BkOQSCeRTqBAazAEiXQS6QQKtAZDkEgnkU6gQGswBIl0EukECrQGQ5BIJ5FOoEBrMASJdBLpBAq0BkOQSCeRTqBAazAEiXQS6QQKtAZDkEgnkU6gQGswBIl0EukECrSGLoag27WuukinW9IJFGgNPQxBv6ut9pBOv6QTKNAaOhiCjtf77yCdjq0ynbF+livQGjoYgo6fONVBOh1bYzqjvRtSoDXMPwQ9P/N0/nR6tsJ0xns3pEBrmH8IFOhSrTCd8d4NKdAaZh+Cy7+pCnRR1pfOiD/MFWgN8w+BI9ClWl86CvTVbW/+Rb/++OO1L6G5y79/H/27qEBp5PxXVX8yrxF/mK+qQGfceO86eBvmOtCFWl863sK/um0Fel0PQ+BOpGVaYTrOwr+2bQV6XRdD0Gt/9pFOt1aYjutAX9u2Ar1uhUMwIukka0zHnUivbFuBXrfGIRiPdJJVpuNe+F+3rUCvW+UQjEY6iXQCBVqDITh69chDOol0AgVagyE4eP2zL+kk0gkUaA2GYO/K2VfpJNIJFGgNhmDvyvV/0kmkEyjQGgzBzrU7UKSTSCdQoDUYgh0FOoR0goUU6MPmy9vbVqDXGYLHsCapdBLpBMso0J9fFejHGIIdR6BDSCdYRIFu+1OBfowh2FGgQ0gnWEKB/vi8aVWg3S538VGGYM9Z+AGkEyygQO82m3/7zzYF2u+Cax9lCPZcBzqAdIL+C3R7/Pnl7+9NCrTjJX8/yhAcuBPpdqtMZ6mLiezejv/+1y1/5I9//+uxWYG+/GI1VjkEQ7gX/mZrTGd5y9n9/Lrrzf3HmZvf/rzxj71SoJsXbvxTX+r5uZEftcYhGI90khWmM+aCyqOW0NX//GGz+cc/d0W4t/3yJgr0o1Y4BCOSTrLCdMZ8pEeTAt0deW6PO7e/bLvz/j3n1J9p8hb+6mXWa7DCIRiRdJL1pbO8h8ptO/PT4ZddEd7tv7lBm89AHYFWJZ1kfeksrkC3Bbj7BHT7y/7jz4db38Mr0I9a3xCMSTrJ6tIZ871mmwL9+XXz7fDLvjk7LVBn4auSTrK+dBZ3BHos0IfDG/mOC9R1oCVJJ1lfOkst0LvN/pfH+xsvBW1VoO5EKko6yQrTGfMs/Ef/X57/ca//9qEATx+Bbuu0y5NIj+6FL0o6yQrTGfM60I//31z+cVd+/25XnfeHm5B2lzR9u+2PbVagq7XCIRiRdJI1prO0O5EOdyDti3P7Pv7Wd/AK9MPWOATjkU6yynSWdi/8/b4/d8V5f/uNSAr0w1Y5BKORTiKdoNm98A+b8yn4G9+/v3vbCvQ6Q5BIJ5FO0P9ydu/ftgK9zhAcWY3pZtIJFGgNhuDAeqC3k06gQGswBHtWpB9AOsHUBbp7Htwrbj6N9J5tK9DrDMGeZyINIJ1AgdZgCHY8lXMI6QQKtAZDsKNAh5BO4DPQGgzBY1jGTDqJdAIFWoMh2HEEOoR0AgVagyHYUaBDSCeYp0B//ofPQNvqYgjmX+vKWfgBpBM0vZXTSaTZ9DAEHay26jrQAaQTtCrQO2fhZ9XBEHSx3r87kW4nnaDxcnYnU6wnokCDDoagjydOuRf+ZtIJmi2ovFuK6X6/HtPdNP2pQJP5h6DnZ57On07PpBM0fazx8WFy95O8g1egyfxDoECXSjpBs4fK7RZE/vF5/0ykbZ16C9/Y7EMw5rO4Rzd7Ol2TTtD6ufDHp3Le+FC5921bgV43/xA4Al0q6QRNC/T0ZI4fn299KNK7tq1Ar5t/CBToUkknaFqgj3eHQ88fn13G1FgHQ9DHWfhXdZBOx6QTtHwu/O69+/7Q88F1oK11MARdXAf6ug7S6Zh0gkaXMR1PvD9s9meR7m9+rvG7tq1Ar+thCDq4E+mKHtLpl3SCdhfS76pz+8u2Oh82TiK11sUQ9NqffaTTLekELW/l/PJ0R6fLmBozBIl0EukErQr07+/71jwuUD/FAagCTQxBIp1EOkG75ewe9h9/ng5FJ6BAA0OQSCeRTmBB5RoMQSKdRDqBAq3BECTSSaQTKNAaDEEinUQ6QbM7kSyoPCtDkEgnkU6gQGswBIl0EukECrQGQ5BIJ5FO0P4z0OPiyuNToIEhSKSTSCeY4yTSNLfCK9DEECTSSaQTzFGgVqRvzxAk0kmkE8xyGZPVmJozBIl0EukEMxWok0iNGYJEOol0AgVagyFIpJNIJ5ijQH9+9Ra+NUOQSCeRTjBDge4uCrWgcmOGIJFOIp1gngvpDyvbjUyBBoYgkU4inWCeAp1kQVAFGhiCRDqJdII5CnSKM0iPCjQyBIl0EukElrOrwRAk0kmkEyjQGgxBIp1EOoECrcEQJNJJpBMo0BoMQSKdRDqBAq3BECTSSaQTTF2gvyylbEHlWRiCRDqJdAIFWoMhSKSTSCdQoDUYgkQ6iXSCVp+BPpzv39xV6hTrKSvQxBAk0kmkEzQq0B+fL27ffHAvfHOGIJFOIp2gUYHePVt/6c5qTK0ZgkQ6iXSCZvfCX75rf/AZaGuGIJFOIp1AgdZgCBLpJNIJ2hTo399fvIW3In1jhiCRTiKdoNlnoBeHoPdWpG/OECTSSaQTtDsLf3rbvvvSWfjWDEEinUQ6QavrQO+fXUfvOtDWDEEinUQ6QbPFRB4mfiLSOLvyx8HH/6DeGIJEOol0goarMR3v6pzogR6j7Moff6y1QQ1BIp1EOoHl7C6cmnOFDWoIEukk0gkU6IVzbyrQWqSTSCdQoE8uanN1DWoIEukk0glaLGf37ddF7bq8E0mBViWdRDqBAj27LE0FWop0EukECvSJI9CqpJNIJ/AZ6BMFWpV0EukECvSCs/BFSSeRTqBAL7gOtCjpJNIJWhfobimRKdaye3QnUmQIEukk0gmaFejPr7ve3PXnVDfDuxc+MASJdBLpBA2fyvmPf+4WVp7ufnirMQWGIJFOIp2g4Xqg2+PO7S/b7rzfXDyic8RtK9DrDEEinUQ6QaMCPS5Cf6xOT+VszhAk0kmkEzR7JtLuE9DtL/uPPz1UrjlDkEgnkU7Q9Kmc21/2zalAmzMER6+eJpROIp2gaYE+HJ8mt+gCXeaJekNw8PqFatJJpBM0LdDTsznvF/xY44VeKmoI9q7cKiGdRDpBs89Avzx9BLqt08WeRGp5s9KYTW0I9q7crCudRDpBs+fCb6vz/nAT0u6Spikey9moQF9+MeGmRmxQQ7BzbbkY6STSCRpeB3p8nvH2ffw0N3O2KNCGCzaNe6xrCHYU6BDSCdo+F35XnPdTPZhzdQU65oYMwWNYMls6iXSCps+FP56Cn+L9+2OTAm24aP3IVW0IdhyBDiGdwHJ2N2l3BKpAJ6BAh5BOoEBv0qxAxz7WNQR7zsIPIJ2gYYHuHoy0Oxc/0XKgKzsL7wh0Cq4DHUA6QcP1QI8rgd5PcxHT2q4DVaCTcCfS7aQTtCrQ43M5twV6t5moQVd2J5Kz8JNwL/zNpBM0KtDdSsq///Xj8+5OpLtlL6jc6k5O14G2I51EOkGjAn3YLwR6KFALKr+PO5GakU4inaDZrZyHmzj3BXpcHXRsaytQ98I3I51EOkHLxUROBbo9BF3wcnbLZAgS6STSCZouZ3cq0EWvB7pMhiCRTiKdQIHWYAgS6STSCWZ5C3/nM9DWDEEinUQ6QbOTSLtjzmOBLnlB5aUyBIl0EukE7S5j+nQq0N01oYtdUHmpDEEinUQ6Qas7kfZXz+8L9H7JCyovliFIpJNIJ2h8K+fBchdUXixDkEgnkU7QbDGR3Rv3SftTgSaGIJFOIp2g4XJ2xwqdqD4VaGQIEukk0gksqFyDIUikk0gnaHUZ02THnRfbVqDXGYJEOol0gmYX0k+1Dv3FthXodYYgkU4inaDZrZxTXDr/YttrK1CrMbUinUQ6gSPQblkPtBnpJNIJ2t2JNNHT4C+2va4CtSJ9O9JJpBO0Ogv/4/Pkx6CrK9CXX3yEIUikk0gnaPYZ6HOWs3uLp3I2JJ1EOoEC7ZQCbUg6iXQCBdqny9JUoFOTTiKdwJ1InXIE2pB0EukECrRTCrQh6STSCRRor5yFb0c6iXQCBdor14G2I51EOoEC7ZY7kZqRTiKdQIH2y73wrUgnkU6gQGswBIl0EukECrQGQ5BIJ5FOoEBrMASJdBLpBAq0BkOQSCeRTqBAazAEiXQS6QQKtAZDkEgnkU6gQGswBIl0EukECrQGQ5BIJ5FOMHWB/rKQneXsZmEIEukk0gkUaA2GIJFOIp1AgdZgCBLpJNIJWn0G+rDZ/Pbn/qtdpU7yhE4FGhiCRDqJdIJGBfrj8+bL+ZuHU5eOS4EGhiCRTiKdoFGB3m0+Xf1utG0r0OsMQSKdRDpBs4fKXb5rf/AZaGuGIJFOIp1AgdZgCBLpJNIJ2hTo399fvIX//a9RN3vYtgK9zhAk0kmkEzT7DPTiEPR+4zPQ1gxBIp1EOkG7s/Cnt+27L52Fb80QJNJJpBO0ug70/tl19K4Dbc0QJNJJVpnOWA8ca7aYyMNTfU5y/KlAo1UOwWikk6wxndEeedtwNabjXZ1TnIA/bFuBXrfGIRiPdJIVpnNqzo83qOXsaljhEIxIOskK0zn3pgK92LYCvW6FQzAi6STrS+eiNj/aoI3fwv/25+P9FNeA7retQK9b3xCMSTrJ+tJZYoEePgHdFehEJ+EVaLK+IRiTdJLVpXNZmksp0OMZpG2B3rmMaQarG4JRSSdZXzrLOwL9+/tm8/tfPz7vrmC6u/VM/P4KqLfvXVKgwfqGYEzSSdaXzvIKdNuBX3Y3Ie0vAb3fXCwO+ra7dy5ir0CD9Q3BmKSTrDCdxZ2FPywfcizQ7eHoDSeS7o/d++Z/o0CDFQ7BiKSTrDCdpV0Huq3M3THnsUC3nfj+9/A/vx4OV7cN+sYnp5vNvwDedr4T6aN/UNP1QE8Fest6oOd/981l7BUo8D4j9Wf/BXpeO/TNw1Zv4YMVvg0bkXQS6QSzvIW/YUHlp6WY33wUnQINDMHRq2tISCeRTtDsJNLu8PFYoNvj0XcvqHys3sen9r3cGkBPBvXmU6Vd+f2H/YWchw7cXRP67ivpT+eQFCjQv2HFea60a/9gf/X8vgN3Syu//yqmywJ9o3abvIUf7/qHtrwN27ty/Z90EukEjW/lPLjhRqR4BPpi220K9OUXy2AIdq7dgSKdRDpBs8VEdm/cb+/P/Bnoi203KNAR7wFryxDsKNAhpBM0XM7uWKG33Qd/eRa+g8uYFOiSXV2FRzqJdILuF1S+uA70jU9OGxTomOtgtWUIdhyBDiGdoPsCveVOpNE3/gtHoIumQIeQTtB9gd5yL/zoG/+FAl02Z+EHkE7Q7FbOywXsbllMpLPVmJyF/4CxnsX9of+DZ78edZFOt6QTNCvQyzfgNxVoX+uBug50uNGexT36/0MP6fRrlemM9VexXYFeHEHeVqB9rUjfQwsM0MEQ9PGzx73wN1tjOqONccMCfVoM5MYCfe+22ywmssj+7GEIOv70o4N0OrbCdJa2oPLPr7/9n91VoMfTQIsu0GWafwh6Pv82fzo9W2E6S3ukx3490LvzG3EF2tz8Q6BAl2p96SzuoXKHBZUfTh+EKtDmZh+Cru9BmD2drq0vnYUW6O5apH11KtDm5h8CR6BLtbp0xvxh3rRADyeTvt30SI8btq1Ar5t/CBToUq0vnaUegR4XFPmiQJvrYAichV+o9aWz3ALd31h085JM79y2Ar2ugyHo4zrQV3WQTsdWmM4iz8IfPSjQGfQwBP3eg9BDOv1aYTrLuw70YiWQ06mksSnQoIsh6LU/+0inW2tMZ2F3Ij3393cF2toah2A80klWmc6y7oVvQoEGqxyC0UgnkU6gQGswBIl0EukEUxfo7srPb8+fyekk0gwMQSKdRDqBAq3BECTSSaQTKNAaDEEinUQ6gc9AazAEiXQS6QQKtAZDkEgnkU6gQGswBIl0EukECrQGQ5BIJ5FOoEBrMASJdBLpBC3Owr+i17Pw3d6t/VGGIJFOIp1AgV7qd72gj+piCLoNt4t0uiWdQIFe6HjFyo/qYQj6/fHUQzr9kk7Q6jPQh/Nj4Q+X1k9gjAJ9+cVqdDAEHf946iCdjkknaFSgPz5vvpy/eTh16bg+vCs9P7XnozoYgo5/PHWQTsekEzQq0LvjE+Ff+260bSvQ6+Yfgp7TnT+dnkknmGFF+sdOn8rZ9ZPLP2r+IVCgSyWdQIE+6XnEP2r2Iej6x9Ps6XRNOkGbAv37+4u38L//NepmD9tWoNfNPwQ9pzt/Oj2TTtDsM9CLQ9D7TZ+fgfZ8muOj5h8CBbpU0gnanYU/vW3ffdnnWfieL7T5qA6GoOMfTx2k0zHpBK2uA71/dh19p9eBdnyp90d1MAQd/3jqIJ2OSSdotpjIw1N9TnL86V74qIch6PfHUw/p9Es6QcPVmI53dU5xAv6wbasxXdfFEPTan32k0y3pBJazq8EQJNJJpBO0Ogs/2XHnxbYV6HWGIJFOIp2g2XWgU1z5+WLbCvQ6Q5BIJ5FO0OxOpCmu/HyxbQV6nSFIpJNIJ3AEWoMhSKSTSCdo9Bnow0TXfj7btgK9zhAk0kmkE7Q6C//j8+THoAo0MASJdBLpBM0+A13CIz1WzBAk0kmkEyjQGgxBIp1EOoECrcEQJNJJpBO4E6kGQ5BIJ5FOoEBrMASJdBLpBAq0hi6GwGIiiySdYJ4C/fkfPgNtq4chsJzdMkknmGM9UCeR2utgCCyovFDSCVoV6J2z8LPqYAg80mOhpBM0fCbS1M/0UKDB/EPgoXJLJZ2g3VM5P+0ejPTp8cUTOsfctgK9bv4hUKBLJZ2g6WpMD4e37vfTPNZDgQazD8FlaSrQRZFO0OxOpC+Puzfy++fJbevUW/jG5h8CR6BLJZ2gWYF+e/rl+FZ+bAo0mH8IFOhSSSdoWqDbQ8/dgej2SHSKpe0UaNDBEDgLv1DSCZoW6OPd4dDzx2eXMTXWwRC4DnShpBM0O4m0P/S8P6yq/OA60NZ6GAJ3Ii2TdIJGlzEdT7w/bPZnke4nWZ1egQZdDEGv/dlHOt2STtDuQvpddR4e7PGwcRKpNUOQSCeRTtDyVs4vT3d0uoypMUOQSCeRTtCqQP/+vm/N49L0kzwkXoEGhiCRTiKdoN1ydg/7jz9Ph6ITUKCBITh69XNY6STSCSyoXIMhOHj9SgDpJNIJFGgNhmDvyrWo0kmkEyjQGgzB3pW7oaSTSCdQoP0a87JJQ7Bz7X586STSCaYu0F+eCG9F+vca9cYdQ7CjQIeQTqBAezXureOG4DGsSSqdRDqBAu3VuIsXGYIdR6BDSCdo9Rnow+ZwG/yhUj3S400jL59pCHYU6BDSCdrdC/909fzDqUvHNcau9LPchQKdgrPwA0gnaPZQuU9Xvxtt2x/flX4WXBv7EUKGYM91oANIJ2i7oPJRr+uB9rTkryPQSbgT6XbSCRTohZ4eOqFAp+Fe+JtJJ2i2Iv3zt/BdLqjc12PPnIVvRzqJdIJmn4FeHILed7qgcm8F6jrQVqSTSCdodxb+9LZ992WXZ+HHPm/zUe5EakY6iXSCVteB3j+7jr7P60D7OgJ1L3w70kmkEzRbTOThqT4nOf5cYYGOyRAk0kmkEzRcjel4V+cUJ+AP217VWfiRGYJEOol0AsvZXejpOtCRGYJEOol0AgV6qZ87kcZmCBLpJNIJFOgza+1PQxBJJ5FO0Oddzt8AABc5SURBVKxAd1cvWc5uNoYgkU4inWCey5gUaGuGIJFOIp2gUYE+bBTorAxBIp1EOkG7Wzknu3zpvG0Fep0hSKSTSCdotpjI5P2pQBNDkEgnkU7QbDm7L6//kzG3rUCvMwSJdBLpBAq0BkOQSCeRTtDsLbwCnZUhSKSTSCdodBLpfpIllF9sW4FeZwgS6STSCRoVaIv38Ao0MASJdBLpBK0upP/xefJjUAUaGIJEOol0gmYnkVxIPytDkEgnkU6gQGswBIl0EukECrSGLoag27WuukinW9IJLGdXQw9D0O9qqz2k069VpjPWX0UFWkMHQ9Dxev8dpNOxNaYz2g9zBVpDB0PQ8ROnOkinYytMZ7wf5vMU6M//8BloW/MPQc/PPJ0/nZ6tMJ3xfpjP8VhjJ5Ham38IFOhSrS+dEf8utirQO2fhZzX7EFz+RVWgi7K+dJZXoM+fiLTZfBt1q8dtK9Dr5h8CR6BLtbp0xvxh3m5F+k+7JUU+Hb6eoj8VaDL/ECjQpVpfOos7Av37+/5O+IfDW/f7aZanV6BBB0PgLPxCrS+dxRXocTGmH59/+/NxX6fewjfWwRC4DnShVpjO0s7Cbwv029Mvx7fyY1OgQQ9D4E6kZVphOku7DvTYnKeF6X98nmJpOwUadDEEvfZnH+l0a43pLOxOpNOh593h0PPHZ5cxNbbGIRiPdJJVprOse+FPh57HJ3s8uA60tVUOwWikk0gnaPZMpH1lPmz2Z5GmeUKSAg0MQSKdRDpBuwvpd9V5eLDHw8ZJpNYMQSKdRDpBy1s5vzzd0ekypsYMQSKdRDpBqwL9+/u+NY9L009xAKpAE0OQSCeRTtBuObuH/cefp0PRCSjQwBAk0kmkE1hQuQZDkEgnkU6gQGswBIl0EukErS6kP7x/P7pzGVNrhiCRTiKdYKYCdSF9Y4YgkU4inWCOAv35tdsC7fZu7Y8yBIl0EukEkxfo/eY1nb6F73e9oI8yBIl0EukEkxfo8crPF6a4kOnju9LxipUfZQgS6STSCaZ/C//QqD9HKdCXX6yGIUikk0gnmOUk0jQ+vCs9P7XnowxBIp1EOoECfaJAq5JOIp3AhfRnXT+5/KMMQSKdRDqBAn3iCLQq6STSCdoW6P6U/CRLMT0q0MgQJNJJpBM0KdDTrZunE/JTrAbqLHxkCBLpJNIJGhTorjb3dx79+LyZskFdBxoYgkQ6iXSCRteB7gp0t6by9tfdL1PciOROpMQQJNJJpBM0uRPpcOf78blIp8XpR+de+MAQJNJJpBNMXqBPT5C7P321bVIPlWvMECTSSaQTTF6gd5vjNfTbA8+nrzpdTGS9DEEinUQ6wdQF+lSW2+PO0yJ21gNtzhAk0kmkE0xdoD+/nt6uXzwNXoE2ZwgS6STSCdoV6N3TqSMF2pwhSKSTSCdoUKCHt/BPH4E+/d64FGhgCBLpJNIJGnwGer6I6dSazsK3ZwgS6STSCVo80mP/zv3uaRnlu15XpF8xQ5BIJ5FOMHmB7u7f/LY/hXR6B3+/meQj0PUV6JhX9RuCRDqJdILpb+W8e/4cj/2CTJ0+0qMvo95XaggS6STSCRosJnJ38SDOw3oi0yxot7ICHXdlE0OQSCeRTtBiObt9a357+nJAfz6846B1dQX68ouPMASJdBLpBI1XpB94An77vr9agY68urMhSKSTSCdYxCM9dp+bKtAPMQSJdBLpBEso0P37/mIFOvYT7gxBIp1EOsECCvRus/m3/6xWoI5AW5JOIp2g/wLdHn9++fu7Av0YQ5BIJ5FOsIAC/fe/HgsWqLPwDUknkU7Qf4HuXCnQzQvTbHwmrgNtRzqJdIKXHfTBElKg43EnUjPSSaQTLLpAX2x7ZQXqXvhmpJNIJ1j0W/gX215bgY7JECTSSaQTdFqgh3vmz2vYK9CPMgSJdBLpBAq0BkOQSCeRTtBpgb6gQD/KECTSSaQTKNAaDEEinUQ6gQKtwRAk0kmkEyjQGgxBIp1EOoECrcEQJNJJpBMo0BoMQSKdRDrBMgr0fdtWoNcZgkQ6iXQCBVqDIUikk0gnWFWB/gugIQUKMNCqCnTGjffO27BEOol0AgVagyFIpJNIJ1CgNRiCRDqJdAIFWoMhSKSTSCdQoDUYgkQ6iXQCBVqDIUikk0gnUKA1GIJEOol0AgVagyFIpJNIJ1CgNRiCRDqJdAIFWoMhSKSTSCdQoDUYgkQ6iXQCBVqDIUikk0gnUKA1GIJEOol0AgVagyFIpJNIJ1CgNRiCRDqJdAIFWoMhSKSTSCdQoDUYgkQ6iXQCBVqDIUikk0gnUKA1GIJEOol0AgVagyFIpJNIJ1CgNRiCRDqJdAIFWoMhSKSTSCdQoLf646DJtsZjCBLpJNIJFOiN/vhjkQ1qCBLpJNIJFOhtTs25tAY1BIl0EukECvQ2595UoCsinUQ6gQK9yUVtLqtBDUEinUQ6gQK9iQJdJekk0gkU6C0uS1OBrod0EukECvQmjkBXSTqJdAIFehMFukrSSaQTKNDbOAu/RtJJpBMo0Nu4DnSNpJNIJ1CgN3In0gpJJ5FOoEBvtcj+NASRdBLpBAq0BkOQSCeRTqBAazAEiXQS6QQKtAZDkEgnkU6gQGswBIl0EukECrQGQ5BIJ5FOoEBrMASJdBLpBAq0BkOQSCeRTqBAazAEiXQS6QQKtAZDkEgnkU6gQGswBIl0EukECrQGQ5BIJ5FOoEBrMASJdBLpBAq0BkOQSCeRTqBAazAEiXQS6QQKtAZDkEgnkU6gQGswBIl0EukECrQGQ5BIJ5FOoEBrMASJdBLpBAq0BkOQSCeRTqBAazAEiXQS6QQKtIYuhqDb5/F1kU63pBMo0Bp6GIJ+nwjdQzr9kk6gQGvoYAhOzdlhg3aQTsekEyjQGjoYgnNvKtBlkU6gQGuYfwguarO7Bp0/nZ5JJ1CgNcw/BAp0qaQTKNAaZh+Cy9JUoIsinUCB1jD/EDgCXSrpBAq0hvmHQIEulXQCBVpDB0PgLPxCSSdQoDV0MASuA10o6QQKtIYehsCdSMsknUCB1tDFEPTan32k0y3pBAq0BkOQSCeRTqBAazAEiXQS6QQKtAZDkEgnkU6gQGswBEevfg4rnUQ6gQKtwRAcvH4lgHQS6QQKtAZDsHflWlTpJNIJFGgNhmDvyt1Q0kmkEyjQGgzBzrX78aWTSCdQoDUYgh0FOoR0AgVagyF4DGuSSieRTqBAazAEO45Ah5BOoEBrMAQ7CnQI6QQKtAZDsOcs/ADSCRRoDYZgz3WgA0gnUKA1GIIDdyLdTjqBAq3BEBy5F/5m0gkUaA2GIJFOIp1AgdZgCBLpJNIJFGgNhiCRTiKdQIHWYAgS6STSCRRoDYYgkU4inUCB1mAIEukk0gkUaA2GIJFOIp1AgdZgCBLpJNIJFGgNhiCRTiKdQIHWYAgS6STSCRRoDYYgkU4inUCB1mAIEukk0gkUaA2GIJFOIp1AgdZgCBLpJNIJFOitXl0PrXuGIJFOIp1Agd7o9RV5u2cIEukk0gkU6G2uPBOie4YgkU4inUCB3ubKU8m6ZwgS6STSCRToTa49F7d7hiCRTiKdQIHeRIGuknQS6QQK9BaXpalA10M6iXQCBXoTR6CrJJ1EOoECvYkCXSXpJNIJFOhtnIVfI+kk0gkU6G1cB7pG0kmkEyjQG7kTaYWkk0gnUKC3WmR/GoJIOol0AgVagyFIpJNIJ1CgNRiCRDqJdAIFWoMhSKSTSCdQoDUYgkQ6iXQCBVqDIUikk0gnUKA1GIJEOol0AgVagyFIpJNIJ1CgNRiCRDqJdAIFWoMhSKSTSCdQoDUYgkQ6iXQCBVqDIUikk0gnUKA1GIJEOol0AgVagyFIpJNIJ1CgNRiCRDqJdAIFWoMhSKSTSCdQoDUYgkQ6iXQCBVqDIUikk0gnUKA1GIJEOol0AgVagyFIpJNIJ1CgzyzzgUfvYAgS6STSCRTopYU+cvMdDEEinUQ6gQK9sNSHvr+DIUikk0gnUKAXzr2pQGuRTiKdQIE+uajN1TWoIUikk0gnUKBPFGhV0kmkEyjQs8vSVKClSCeRTqBAnzgCrUo6iXQCBfpEgVYlnUQ6gQK94Cx8UdJJpBMo0AuuAy1KOol0AgV6yZ1INUknkU6whAK932x9envb7oW/zhAk0kmkE/RfoD+/bvZ++/OtbVuN6TpDkEgnkU7QfYH+/X1fnT8+b37/641tK9DrDEEinUQ6QfcF+nA89Nw26Lc3tq1ArzMEiXQS6QTdF+jd8chzeyT6xsegCjQwBIl0EukECyjQTy++uLptBXqdIUikk0gn6L5AT35+3Xx5Y9sK9DpDkEgnkU6wmAJ9eOU0/OaFyTa+fIYgkU4ineBlB32whCbrsO0B6K9n4RXouxmCRDqJdIKFFOi2P//xz7e2rUCvMwSJdBLpBMt4C7/tzzevo1egiSFIpJNIJ+i0QH98PhwPHy793N2M9MZFoI8KNDIEiXQS6QRLKNDtN28ffyrQyBAk0kmkE3RaoJfe2Z8KNDEEiXQS6QT9F+i2P988f3TYtgK9zhAk0kmkE3RfoO85/37ctgK9zhAk0kmkE3RfoPfve//+qEAjQ5BIJ5FO0HuBnlYDfceaygo0MASJdBLpBL0X6Ol0vAL9GEOQSCeRTtB7gd6ybQV6nSFIpJNIJ1CgNRiCRDqJdAIFWoMhSKSTSCdQoDUYgkQ6iXQCBVqDIUikk0gnUKA1GIJEOol0AgVagyFIpJNIJ1CgNRiCRDqJdAIFWoMhSKSTSCdQoDUYgkQ6iXQCBVqDIUikk0gnUKA1GIJEOol0AgVagyFIpJNIJ1CgNRiCRDqJdAIFWoMhSKSTSCdQoDUYgkQ6iXQCBVqDIUikk0gnUKA1GIJEOol0AgVagyFIpJNIJ1hVgQK0NW6Jjfqn3bhtgMbGLbFR/7SulP2EwI6XU3bPZ9/xubc/odmznYsdL6fsns++43Nvf0KzZzsXO15O2T2ffcfn3v6EZs92Lna8nLJ7PvuOz739Cc2e7VzseDll93z2HZ97+xOaPdu52PFyyu757Ds+9/YnNHu2c7Hj5ZTd89l3fO7tT2j2bOdix8spu+ez7/jc25/Q7NnOxY6XU3bPZ9/xubc/odmznYsdL6fsns++43Nvf0KzZzsXO15O2T2ffcfn3v6EZs92Lna8nLJ7PvuOz739Cc2e7VzseDll93z2HZ97+xOaPdu52PFyyu757Ds+9/YnNHu2c7Hj5ZTd89l3fO7tAyyWAgUYSIECDKRAAQZSoAADKVCAgRQowEAKFGAgBQowkAIFGEiBAgykQAEGUqAAAylQgIEUKMBAChRgoFUW6MPmy+GLv79vDn7/a97/oxbudzv66fTdw7PvVu1yx0u94o932x39xz9P3xV6yS93fN6XfI0F+vPrqUC3X1UZp9Ou/vbn/tu7w3dP07Vaz3e80Ct+3tfjX/Z6L/mXZ98p0JHsAj3+nfrxucAc7W1/DO8a5Mfnw1+j+30Ep+/W7MWO13nFd3u+K8vtS/1t922pl/xyx+d9yddXoNu/QucCfSjyjma3p4cjsO3ef3s6Bj98t2rPd7zQK37e87v9Lhd8yQ87PvNLvroC3b6R+bf/PBXo/emL1bs7Hnhsfzp/2v2dOr6Ru1t9nzzf8UKv+Pm1vd8nUOolv9zxmV/ytRXo9gfwl+0snT4WOn4iuH7nsdl/cWqV7V+utX8k9nzHC73iZ4ceKfSSnxwLdN6XfHUF+u9/PZ4L9OfXf/zfr2XOTO7t38kdj8Yen97trN/hLWzBV/zwYhd8yY+7PPNLvrYC3TkX6P7j0CJnJo/24/N0CP7jc5FpOvZGvVf8eEal3kt+3PG5X/JVF+jD8dKGu/WfmTzaHodtd/XpOq4q03Tc8XKv+MPp6p1qL/l5x+d+yVddoKePg7a/sfYzkwfbMdrt8eU01djz446Xe8UfTlfAVnvJzzs+90u+6gI9W/+Zyb3tFO3/TlU7HDnv+JMir/jO3W7fy73kj8cdv/x2jpe8RoHel3hHt7uDYP9DuNoHYucdf1LjFd/bn0yp9pLvPJ0425vnJVegq7F973YcnctTsgXOpjzt+JMSr/jRbl+LveQHz19kBTqWmgV6WSMXFwXW2vGzCjt+st/XUi/5kQKdxqlALy6oX//9KdsauVyXp8xtKc92vNIr/rSv97VuPnu+43O/5Gsu0PPPpJ9fV/+x0Ok09Pm7KjdGP9/xQq/4+fTz6RaCMi/58x2f+SVfdYFu/zJ9evzlw+ZVun/+NrbO0jzPd7zQK77rj+O+FluA6/mOz/ySr7pAD1eLVbiz77wm4mlnqywO+XLHy7zij0/34Jzfu9d4yV/u+Lwv+boL9Lha9drf01zcz3b+i1RkefJfdrzKK75392xfi7zkO892fNaXfI0FCtCEAgUYSIECDKRAAQZSoAADKVCAgRQowEAKFGAgBQowkAIFGEiBAgykQAEGUqAAAylQgIEUKMBAChRgIAUKMJACBRhIgQIMpEABBlKgAAMpUICBFCjAQAoUYCAFCjCQAgUYSIECDKRAAQZSoAADKVCAgRQowEAKFGAgBcqkHjbPfXl8/Pv75h//3P/Dn//r+G+dv7jm59fN73/dsNmnbcCEFCiTSgV6t/vu2RdXKVC6pECZVCjQh82xN89fXKdA6ZICZVKhHBUoi6dAmZQCZc0UKJNSoKyZAmVSr5TjodxOH45+O39x+vfPX+9sm3Oz+fSiQLffnevx/PXd/r/87c+Lbex+Of7G44/P5z/h5TYO338aca+pQoEyqdsKdNtyB6d+vD98+/t/PT8CvX/qv+MG7s+nqfb/XijQl9s4f//mYTC8pECZ1E0F+tRlx3Y71+J/e16g23/xdMB4t6/Ip/48bO96gb7cxv4Q99lBMLybAmVSVwv0lc9At//g0HcPx3fUu3Lbtd7D+cjy6Y84FuPhvf3uX/zyePFfXC3QX7ZxfyzS3T+46WNWUKBM7Pl1oPvDxqsF+nCusG0PHo8rn952P2u3++f/7cP5I8ztv3lsztcL9Jdt3J0OPJ134nYKlEndUqB357Y7/NZF/d2/KNDzSaWL/+b0D2KBvtzGL38y3ECBMqkbCvTyTPv+Q86Lk+0X59APjkW4/Vcuz57v3ofHAv1lG6f/QR9/MoQCZVI3fAZ6eTpn/5nnRdv9ch3o08n3b+ffuDg3FAr0+TZOlz85Cc8QCpRJ3VCgF+fHD+V2cdj5S4EefuP8ueXlf5sK9JdtPF6ewVeh3EiBMqnbCvS1knz55dH+PfzpcqbTgeU//vnGZ6C/bOPx8j/XoNxIgTKp297CPz8Nnj4D3f3Ol/M7+PvzZfFvFOgv27jcmNNJ3EqBMqkbCnR3Bujbi3/x9Bu/nivf/sPf/9/3831Hp558ePYW/uKk++FPeG0bp0a9Xq5whQJlUrdcSP90HPnyEqPT9fGX7je//e/P5yvwn66rf1Ggx7P0pz/hl23cnf8HFSg3U6BMKhfooRRPX+zfRe+PD+8Otwud6nB/6udlge5/89ibd8da3J8P2v3eaRv3xw82z/cy/bKN0x1J+39iQRFuo0CZVCjQwxnxbxdfXF40ero36eC/v7Kc3d1Tqz672vSyQJ8uW/qfxz/hl23cPX3vAJQbKVAmFQp0f9X78Rkfm5dXMn05/+f7Y8f/eqVAHy5Om59b8H8cPuQ8f7R56ssv5/P4L7dx2PqpeeEWCpRJhQI9VNenyy9OF2V+uvyXf10P9ODZZ5aHWjzfFn9xbujuUI0Xf8LLbRwPU92MxM0UKMBAChRgIAUKMJACBRhIgQIMpEABBlKgAAMpUICBFCjAQAoUYCAFCjCQAgUYSIECDKRAAQZSoAADKVCAgRQowEAKFGAgBQowkAIFGEiBAgykQAEGUqAAAylQgIEUKMBAChRgIAUKMJACBRhIgQIMpEABBlKgAAMpUICB/j/J3izY5LR6EgAAAABJRU5ErkJggg==)
qqnorm(fit3$resid) #check normality of \epsilon
qqline(fit3$resid)
#refitting model using ML
fit4 <- lme(nitrogen~1,random=~1|influent.f,method="ML")
summary(fit4)
## Linear mixed-effects model fit by maximum likelihood
## Data: NULL
## AIC BIC logLik
## 262.557 267.3898 -128.2785
##
## Random effects:
## Formula: ~1 | influent.f
## (Intercept) Residual
## StdDev: 7.159246 6.534319
##
## Fixed effects: nitrogen ~ 1
## Value Std.Error DF t-value p-value
## (Intercept) 21.2171 3.165066 31 6.703524 0
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -1.93438234 -0.55703907 -0.03545401 0.83759580 1.93232419
##
## Number of Observations: 37
## Number of Groups: 6
#results should be similar to fit 3
##using another package lme4, lmerTest
library(lme4)
## Loading required package: Matrix
##
## Attaching package: 'lme4'
## The following object is masked from 'package:nlme':
##
## lmList
library(lmerTest)
##
## Attaching package: 'lmerTest'
## The following object is masked from 'package:lme4':
##
## lmer
## The following object is masked from 'package:stats':
##
## step
fit5<-lmer(nitrogen~(1|influent.f), REML=TRUE,data=ex.data)
summary(fit5)
## summary from lme4 is returned
## some computational error has occurred in lmerTest
## Linear mixed model fit by REML ['lmerMod']
## Formula: nitrogen ~ (1 | influent.f)
## Data: ex.data
##
## REML criterion at convergence: 252.4
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.91444 -0.53646 -0.03217 0.83714 1.95824
##
## Random effects:
## Groups Name Variance Std.Dev.
## influent.f (Intercept) 63.32 7.958
## Residual 42.66 6.531
## Number of obs: 37, groups: influent.f, 6
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 21.223 3.429 6.189
anova(fit5,df="Kenward-Roger")
## anova from lme4 is returned
## some computational error has occurred in lmerTest
## Analysis of Variance Table
## Df Sum Sq Mean Sq F value
![](data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAABUAAAAPACAMAAADDuCPrAAAAZlBMVEUAAAAAADoAAGYAOpAAZmYAZpAAZrY6AAA6Ojo6kNtmAABmAGZmZgBmZmZmtrZmtttmtv+QOgCQtpCQ2/+2ZgC2Zma2tma2tra2/9u2///bkDrbtmbb////tmb/25D//7b//9v///+3kuRjAAAACXBIWXMAAB2HAAAdhwGP5fFlAAAgAElEQVR4nO3daWPiVsCeYZK04868XTwNXUKDPfz/P1kQi1YwPCMBOlzXh3hsYxBCunO0IBYbACKLR08AwFwJKEBIQAFCAgoQElCAkIAChAQUICSgACEBBQgJKEBIQAFCAgoQElCAkIAChAQUICSgACEBBQgJKEBIQAFCAgoQElCAkIAChAQUICSgACEBBQgJKEBIQAFCAgoQElCAkIAChAQUICSgACEBBQgJKEBIQAFCAgoQElCAkIAChAQUICSgACEBBQgJKEBIQAFCAgoQElCAkIAChAQUICSgACEBBQgJKEBIQAFCAgoQElCAkIAChF4yoL9+LhaLP/9pfPNt/Af5eNve7/vAL5aLyqWHPH+Tzx+L2vEpnH2oo9WlXwKp1w3oMU/3DWgzgMf+3XKTVkCP0/1FQLd/I6AwhRcO6B9/19/cK6DrVv4Gu3b5Jp2A7if8ckBXXwxPgdQLB/RQzXsGtBPHY8NvuEk3oNWvLwb0y+17IPXKAd236Y4BPdTv++7fq4ub6GdvUv2++u0htd+HH+ridADjeOWALv76d3PXgK6aQdyn8nvnr766SSOg+0fYPQcBhcd46YBWXWkHdFWXdWd3OPy9itYff6/3v1ge/7D6cT0+PG5cH/+2H67D/Ry/rW7QGYJ+eZNmQKsp7we09RRWp439bqqB3/baAd1VphnQOjeHH+xq+d9/7EeFVUD/3+kvl40IN+7ylKp+QNeNOz7eeWcX55c3+TKgnacgoDCh1w3of/qxr04joHVtjiO4ZaM/3cM7B/u+LZs/qlrWD+hq0f7Jut+1L2/SDOj6OJ2Nh+o+BQGFCb1uQL8dNsnrgK5b1al+tGwMMtfHcV3Vq6pIq2OZqh993zSHkP2AdoeT1S3ae1+/vMkXB5EGnoJ9oDCZFw5o9d/vdUAbQ9H1sWPLxtjtFKf9T3e3rHJ2/JN6zFr9azigzR2ajb+++iaXT2MaegoCCpN54YBWjfnzn1N0TslphKg5Iqw3p1fHnw4cwT8f0NYbSHf6Af36Jt2AVtNzeqihpyCgMJlXDuh+CNqK5fEQ+uoQsuaIcH3aml+dftr8m5391v39RqCd41VDT0FAYTKvHNCqLX/+32ZAj61aNwJ6LFI9vBsMaOM4/IWADu3gPB2d+nblPtD69o1bHQPafQoCCpN56YBWwfnPh29am+PHWF4d0MOBpYsBrQ+xf/5439SHoBoBPXeTk+ZR+JPjQw0+BQGFybx2QE/Z+90R6PGO3i8dRKqP0K9PO193t2gE9NxNTi4G1AgU7uu1A3o6TXJoH+ju39cGtN70vhDQ09uM9l38b4vD3TQCeu4m7fu4GNDuUxBQmMyLB/S4S/HSUfgrAtrI2oWA1m90r/dkdt+E/9VNLgfUUXi4qxcP6HEI2jupc1GfB3pTQKt/nQto41JLh1Fn73p2X93kckCHnoKAwmRePaCHYNXb8/13Il0R0NO+younMfUv9jnw/sovbnI5oENPYfAvgDG8ekDrN2hu2m9o77yvaHN5H2g7e9Xvrrki/dBl9C7f5IuADjyF4+lVGgqje/mAtj4fafBqTFedxnQK1//4eYjZ15+JdKZrF2/yVUD7T+FY5AkueQqv7uUD2rmC3Lodm+tPpD9d/Ph4KtFXn8r517/7x7rtJl8GtPcUjvsVBBRG95IBfRa7kDffCJrdBHgUAX2oj7cvj49fcRPgMQQUICSgACEBBQgJKEBIQAFCAgoQElCAkIAChAQUICSgACEBBQgJKEBIQAFCAgoQElCAkIAChAQUICSgACEBBQgJKEBIQAFCAgoQElCAkIAChAQUICSgACEBBQgJKEBIQAFCAgoQElCAkIAChAQUICSgACEBBQgJKEBIQAFCAgoQOhPQ1WKx+H7fKQGYmV5Al3/8ve/nYvHtERMEMBedgP76uVi8bz5/LP76d/uf98dME8AsdAK6XuwGoOtdRbf/NQQFOK8T0GU16lwt/vxnsx2H/vXvQ6YJYBbaAf38sSvndjv+2+kbAIYNBfTzR3UEXkABLhkK6H5H6ObjzSY8wHntgP76uUvnstoFuv1y80Gk3UH8PYNXoHidg0irxV//frxV5Vzux6E3WC5aHMMHyjZwHuhiV87PH7e+FWn3p3Vyd39vDwBQtO47kXYZ3B+Jv3UjfN0ZsW4HstcVeEGpbluCYHbGW8aX3V6urxuCPnolZ0KjLV3wlEZbxD9/dHeZXnkalLWsCEO99NJSuv4iftiXub5xD2Y/lwL6OobHnF5aSte/GtNuLdgGdLW47VoiRqCv7FhOAeW1dBfx3clL1dmg24HoTYeRtrdP94He8Cg8p+OLKKC8lt7VmLbV3J9O3z8qdNkqPApvLXtyVxwNqn/fuKWDSJRv6GpMh4De+FbO4wlQ9XfX/bm17KlddUR9IKCOwvMKhq/GVAX05ouJZO9EspY9teG9m70btf91qKeXltKNGdD9+4/2rv9Ta9lTG967eeZG1b/qAauXltJdCOhUV2NytvWMDO7dPH+z9kvqpaV0F/aBria6HIiAzsjVAW1uul/1J1CA80fhP97u86Fy1rJndn1AB/536KWldP3zQBfvVUCXt1+PrrqUU73zc3XdjlBr2TO7MqDDR929tJRu+J1ItxxFP1l3/k5AS3DlQaTBXTFeWko3sIjvE3rzAaTDdZhXp0GogJbg69OYzu/J9tJSujEvZ7cfem5Dug+ngJbgq2N9l3833XTBMxhtET8cut9UJ4MejuILaAnCem68tJRvxOuBnnJ5uA6JgJbt65PQvLSU7rSI128iqt3yTqTm+5b2BRXQkl1zCq+XltKNFtB6E776ZvHXvwJarCvfAOGlpXTjH0SqbHP8TUDLdP37x7y0lG68RfxwGlP93ZUjWGvZrNzy7lsvLaUbcRFft06+3+0SENDC3FLPjZeW8o26iLevYb8W0KJcv+le/8VkEwNPoXkQ6a9/O0eSbr4eaDQF1rI5uLmeGy8t5RNQvpbUc+OlpXyPX8StZU/u9k33+i9Hnxh4Ko9fxK1lTy2u58ZLS/kev4hby57X79Rz46WlfN3PRGpcxW59+yXtoimwlj2nfNO9vofRJgae0sWAOoj0un67nhsvLeW7ENClgL6efTbHqOfGS0v56kV81b+YyDSfytmdAmvZ02i++KPc3Qh3Ak+sXsSrj4RrucsuUGvZExk1n15ayndhE/5eU2AtexrHeAooXEVAOWiOPcd5Uby0lO7xi7i17AEam+mnQWdzy11A4RqPX8StZXfX2NXZOmwkoHCb3iK+dDGR4jU21FuHjQQUbtNZxNuXYxLQMtVb76cdn/sfNX4x5uNAqTqL+GrRuibyXabAWnZnzYFm+7DRuMeQvLQUr72I//p5n5PnW1NgLbuzUyWbx42OAXUeKFyvexrT4v3uU2Atu7NWK1s/bB2eH+mBoGAC+nrqdo5+4lL3gca/T3gmvU34e+8CtZbdW+uo+8iHjbqPNMGdwhPpLOJ3uoRdawqsZffUeb/RyIeNuo81wZ3CE+ku4svFH3/feQqsZXfTO3t+7MNG3Yeb4E7hifT2gToPtFTtseei3nifpp4bLy3lE9AXMdkw8+Jj3vfx4N4ev4hby6bXOOp+z9ntpaV0j1/ErWUTW7Tyec+Gemkp3eMXcWvZpNo7Pg9f7/bYd3ogeJBLi/jnf9gHOm+DbzYSUBhLdxFvHUaa5iDSomuKB2HTmL+dGX2vWe6lpXQXL2c3zds6BfQOejO5OacFFMbRu5zddtS53l2S6ePtThdmspaNr53M9peNgMJYeu+Ffz9+tNznj/u8J8laNrI6nQIK0+qeSL/b7fnrZ7Xzc32fIai1bDSL3sizdckQR+FhZEMB3SyrseedPuPYWjaS3v7OQzkFFKYyGNBVdfjoMBCdfAqsZbHhQ0WLdjmbAb3vYTsvLaXr7gOtxp6r6qqgAvrkzoTzbEDv/n54Ly2l6yziyyqdH2+7jff9f6efAmtZ5LSd3qxnfezo+LvNAzbd60m88wPCnXUW8Y+3XUGrC9Pf6wPmrGU3O7/dfvq4o2NADzc//Nndp/PODwh31l3E959rvK5WyPtcnN5adrVeMTufz7Hp/Kyf2HtP750fEO6st4h/vO0OIO3ekXSXDXhr2dXO7O48G9DGLx9Rz42XlvI9fhG3ll3psJF+PDB06mSdzXrj/VHJ7E7xo6cApvX4RdxadtnAqLMTz964c/Mkc/UpJgIm9PhF3Fp2Xm8nZnOz/VTKp9jfOeQpJgIm1D2RvrUOOg/0oQby2QjophXOTTuij570vWeZDpiKgD6vuov940KDg85nm5XPNj0wtrOL+Mfbfc5ispad1ezk6RhRPfrsdfXR09vzhJMEozq/iB/eFz/9FFjLhrWOtjfi2SnnoyfzgqeeOBjBhUV85Z1ID3UhoHOo58ZLS/kuLOJr+0Afqnu+51Pv7hw0i4mE3yCgT6txllJjv+dc4rkzmwmF0IVFfHmfN3Nay87pRrO+Lt1MzGtq4XbnF/H9ZUXuMAXWsmGz3XI/mdfUwu0unQfqKPwDzb6eGy8t5bsQ0Pt8qrG1bNhsq9kw52mHazx+EbeW7Z0ZcM559sx52uEaj1/ErWU7Z85VWsx6/sx40uEqj1/ErWU7i6E3vG8EFJ7awEd6nI4fff6oPqNz6imwlm2Gz5gXUHhyrUX818/G6vu+Wd7lOJK1bNO5avLmeOXkGZ762TLjSYerNBfxj7fjByEdSnqPAai1rHfN+dPbj/a/evTk5eY87XCNxiL++aNx5md1QpOLidxBZ7NdQGE+Gov4snXm/K6g73eZghdfyy4E9LhJP1dznna4Rr2Ibzfgm8HcHU3yVs7pnbrZPYI093puXv6l5QXUi/iqde2QahPexUQm1jxu1A/o7OfNzCcfvnRaxH/9bA0414tvKwGdVOfA+7BHT+Rvmfnkw5dOi/h2yNncgl8u3l0PdEpnNtnLqefmZV9aXsi5gG6iCyrXJ5Je/6evuJad+lgPP+d+xH1Iac8HusYM6LI9frryJKhXW8sGBprHfZ6lzYrSng90ndsHuukeVfrSbvRZn3m/Owh13Z+/1Fp2bldnKRvtbaU9H+g6dxR+qKiXrTtvXPp4u+7vX2ctGyxmsfXcvNJLy6s6ex5oL4hfWXZ7ub5uCPoqa9mZIWcZZ3wOK/JJQcOldyLd9E7O/qWbPn8M70MtfLu17UwxN/3vHz2lEyjySUHDhffC33YWaD+XAnp5yNmu6aOndQplPiuonb8a042H4K8fgXanoMS17Hw6u9deWggozFZ7EW+ciHTrlZj6x5xedx/oUD4vBLTYYXiZzwpqZ65IH1xGZOUo/NFAOs9F9HjduiIV+8TgYLxFvL3Zv/vuFc8DHRx4tuo5OCwtU7FPDA7GXMRf/p1Iw5vu5+tZ/7JMxT4xOBh3Ea8ugld5wffCXxh7dr95hXpuCnpp4YzHL+LFrGW9Tm56ezxLL2bHCz1VXtTjF/F5r2VnB52b3qiz6ONFg17r2fKKHr+Iz3ctO7tRfsmjJ/qeXuvZ8ooev4jPdS0bGnYeN9zVs/JyT5iX8/hFfK5r2eK0g7Pd0t7nbL6ul37yvITHL+IzXcsGBp8DB40ePZWP9erPn/I9fhGf5VrWbGT9ZVP/aJ5Pa2TmAaXrL+KHK8uv7/KRnJvZrWXndna+0EWmrmYeULreIl69nWgb0NWi+xlJU03BnNaycweHegF99IQ+BbOB0nUX8eU2nr9+bgO6HYje5VON57WWDQ412yPSR0/iEzEzKF1nEa8+ibMK6MBndEw0BXNay05Higa23jd2fXaYGZSus4gvd9vth4B+vN1lN+ic1rKhnZ92fZ5lZlC69iK+v4j8IaBXXlH+t6dgNmvZpWyq5wDzg9IJ6PXs8ryRuUPpLgTUJvxRf6Bpl+c1zB1Kd2Ef6Ormz0XKpuDZ17KBeNpov4q5Q+nOH4X/eLvPiaBPv5Z1Dht5p+bVzCBK1z8PdPFeBXR5+wdzhlPw7GtZ4y2aRp83MYco3fA7kSr36eczr2W9Y0bNnz164mbATKJ0A4v4PqH3eiv8865lvdHm8T9PO8XPxoyidI9fxJ92LWuPPjcbh95vZUZRuscv4k+7lnXOWGoXlSuYUZTu8Yv4061l/S339qWSHz1982FeUbrTIl5/pHvtFd+J1JkFm+ZpTJunm9rnZmZROgHt6Bx43wwlleuYWZTu8Yv4c61lp0zWAXXmfMr8onSPX8SfaS07xXNT7/rc1CXlNuYapXv8Iv40a9nAlvvmmNOnmchZMdco3fCHyi3u9j6kZ1nLBo8d2fX5e8w1StddxA/5vGNCn2Eta5fydLxdPX+PGUfpOov4rp/7azCtXuZiIs0t9+N2+2k/6IOnbd7MPkrXWcRX9XvgP38Uezm7eqTZPWfpdAR+IaC/z+yjdO1F/HAp5b11oRdUXpw18Ms7T1tZzD5K1/1Ij8ZFmEr9TKRF9+zOeuO9vefT6v+bzEFKN/SZSEPfTDgF9w/oprXjsxtQK/5YzEdK192Eb+z3nOpD5Ya2nO+jWc3m2PO0+9NZ86MyIyldZxFf1+9/38b0+zQP+aCA9h62cXnPReNjO6z3YzEjKV13ET+dvbS+10Xp7xnQ7jlLm+bw06GjsZmRlK67D7Q3PJz8kkx3DGhjlHkhoHeanBdgXlK61wlocwfo4fvGF9dKnoAZSukev4hPv5Z1RpfnAjr1ZLwe85TSPX4Rn2gta/ayvW3eH2TbdJ+GeUrpHr+IT7KWNbJYf23/sv2OzgmmAbOV0j1+EZ8qoMevxz6eKmm1vhdzmtI9fhGfKKCHL+1N9+kekAHmNKW7cD3QGX+oXOu6IPUP7e28L3Oa0nUW8Y+3e57BtJ+CSQLa2APa+eH4j8YZZjal6yziy8V9LgLanILJAnr415QPxCXmOKXrnUg/zfvfL03BWGtZ42hR+5ylkR+Ha5nllK4X0HsPQMday4bP6mxsxFub784sp3S9y9nNdQRan9LZPGfJ9ZEfySyndL3L2d3nEkzNKRh7BNr7Tj0fxFyndP3L2d3l0HtzCkZZy5q5rH8wxj0T8wJQujLOA+3s9uQ5eCkoXQHngXb2cFprn4aXgtL1DiLN7DzQ+oT5ce6PEXkpKN28zwNtjD1twT8frwWlm/F5oN0td6d7PhuvBaXrvZVzJiPQ3qmdTvd8Pl4LStc7iHTvs5iitWywlOr5bLwclK67iH+8/fH3nafg1rXMOHMuvEqUrrcP9LdOY6rPIr3+T29by5qTd+PEcW9eIko3ZkCXrT9efLtyCm5Yy07nLGnoHHh9KN14i/hu9Flv/u9SfN3b6q9ey+pkOuI+D14fSjfeIr5etPeefrxdd0T/urWsNeJ0zuc8eH0o3XiLeO8MqCuv7HTNWnbuzZrW0Kfm5aF0lxbxz/+4YR/o54/u4fvPH1ftQ/1yLevt7RTQmfDyULruIt46jHTLQaR+LkcJ6NDBIgGdCS8Ppess4u3D8Le8rXOSEejwoXYBnQkvD6XrLOLV9ZTXuzOQPt6uPQ9pr/9pIL+7D/T8iUoOIs2D14fS9S5n974bOu7K1x9SXrYa6Sh861TPc3/iNKY58PpQuu6J9Lut7l8/q23v9c1D0MZO0913yXmgV73TyIn08+D1oXRDAd0sq7HkfiB6gzHeiXRlG9VzDrxElG4woKvq8NFhIHqL+hhU+l74Yz2tewXwIlK67j7Qauy5qvZeBgGNpqB3Uc/ej5knryGlG7yg8sfbbuN9/9/pp2Bx+sfw242YK68hpet/Kuf3wylJ2//cdBApnoLj9UEEtDReQ0rXXcRXi11B1ze/EWmnuhpo/Uer6+5gt5ad2imgJfEaUrreIv7xtjuAdP3V6GrrztH3cwFdDDr9sv2VGfMiUrrxFvHDW5dWp0HoLQFt/vLwdbQJ41G8iJRuzMvZ7Yee25Duw3nDJnzzO+fIF8OLSOlai/jx3evbBt52JZHN6QyoTbX5fzgZKrqYiHoWw+tI6RqL+Ondl4d3FN22D7Rx7aXtHe3+fXVAKdYtSxDMT2MRXx5GnevqQPytV2NqXrxuX1AB5ZYlCOanXsS3xaz6ua1fNfa88RPi6034431cGdDfMdUaOtma746hJPXqsTpss3/+2F+G7tYT6ZfNm2/v5JuAumMo22n1OF0QeX08frS6bS9oe5u/Og4loO4YSnZaPbZjxmM3D9vi6xsDuDuRvk7o7lR8AXXHULJ+QJfHgeetAe1+svHtf38zAZ3vHUMBegHd7b3c/+QOAfxdAjrfO4YC9PaB1h9ldOM+0EcQ0PneMRSgXj0OR9FPu0B/43J2N54BlRPQ+d4xFKBePfZb7Id3EW3q80IDAuqO4RXUq8f+BPr99UA3zZLeTkDdMbyCxupRnbq5aFxPJB2ACqg7hpfQWj12VxGpxp/Vp2vG/RRQdwwvYXj1+M3PQxJQdwyvYJLVQ0DdMbyCea8eAjrfO4YCzHv1END53jEUwOoBEBJQgJCAAoQEFCAkoAAhAQUICShASEABQgIKEBJQgJCAAoQEFCAkoAAhAQUICShASEABQgIKEBJQgJCAAoQEFCAkoAChuQd0vdh5n+Kuf/38859R7/Dzx3ZSR77Po9EntjLh3IUSzDugv34u9v76d/w7X40cu9VhWifp0dgTuzPp3IUSzDugy8Xi+6ZK07fR73s18mhxXU3kNkp//D3ivR6MPbGVKecuFGHWAf14OwznVuNXaTny5va2nN+aX0c19sRWppy7UIZZB3R13Lj8/FGNlcazjcfiv4zapFOO1qOnbvyJrUw3d6EUsw7o8jiY2w7rRl3Fd4d73sct3Wkn5fa+x90LOsHEViabu1CMWQf0ZOwx0uePb2MPFZfH8dzoOZpgYrsPIKAwqIyArqY4tj1qk+pdn5PsBJ00oJPMXShBEQHdDpEmONNmooDWm8ZjmjCg08xdKEEJAd2drzjBEGnsgB43g2cW0InmLpSggIDu1vAp9tEJ6M5UcxdKMP+A7g5CT3Kit034zYRzF0oww4DuTnvcOpzcvftunDW8fb8OIu2MN3ehRHMP6O6bkbYwJw3ohKcxVaYJ6IhzF0o0w4A2rSe8WNBcTqSvTBLQKeculGDeAd2NkCZbw8dt0sfbYWg7zVhxinuddO5CCWYd0OpNjJPd+7hNmvRiIpMEdNq5CyWYdUBXk+6hG7lJ6+q6mhNdzm6KgE47d6EEcw5odYX3o/GHdWM3adILKo8f0InnLpRgzgE9HDafS0D3SZroXZHTXCRPQOGiOQcU4KEEFCAkoAAhAQUICShASEABQgIKEBJQgJCAAoQEFCAkoAAhAQUICShASEABQgIKEBJQgJCAAoQEFCAkoAAhAQUICShASEABQgIKEBJQgJCAAoQEFCAkoAAhAQUICShASEABQgIKEBJQgJCAAoQEFCAkoAAhAQUICShASEABQgIKEBJQgJCAAoQEFCAkoAAhAQUIvVhAl4uW75uPt8X7GHf8v/7Z/ueKOxu4ya+fh2n5/Yce69kAVxHQMZKzWvyZBvTzR2NqfvehBRTu6sUCurde/PH34Z+jBvQK3cfbfr/4tvvHLqTffvehBRTuSkAfGtBtNk9/uY4KKqDwOAL60IAum3+4qicrfGgBhbsS0MV7dRTnr38PP6l2StYl221jH7/9ePvzn91e1Pf2zdb7XZjvdb/WrV2a+x2v+861E9f+bjsd35o/PE1m9w6qx64muPPQpz9tP4nqOeRHqYAzBHTx3xaNPp0OMx2Cevz22/7Gf/7Xwy2bN+sG9HBY/XAX69NBovdNN5mdkeu6+rYb0N4dHCZ498szAW0/idWi9QyBsQjoPiyrQyNX+5HacX/kcv/18GV34+1v/0/vZu2j8Msqdev9TbZjweNYsVXHncOQ86QzjNxPZv8OOhPcPwrfnrrtD7/vvyQHqYDzBPQwMFsOjf4O7Tnuntx+ux/TdQeJ3Yp935x+eBpknu6iFdDWdvW2ld97d96/g+YEDwW0dwfHgaghKIxLQA8NWx9yd9wXutyN107NOaXt2NPWzdoVWw3u7qyTXP9wf6+1fVD7+0A7d9Cc4KGAnnsSwMgEtFmrxpBwV6HGJvZyF6XjjTs361RsObSpvNvy/q2A9u7gzOD3vTd1aweQYCIC2uxR421Buz2N3RQdb9y5Wbti3R2bh2Pgh6M+t2/Cn7uDswHtTd2yeVQMGI+AXghoY4RYderagLaquDy0b2AT/qqDSOfu4PqAno7jO0cUxiWgnYA243dhBNreJr4wAl3VpzP1Ato/jWkgkefu4FJA+1vsawWF8QloZx9oM34X9oG2t9IH94FWt6pvOnAUvn8i/S6V5yZpdf0+0MHzldY242FkAtrq0bI9JOwfhd/fuHOz4aPw1d+choPbrPUDerqj6qTRw+mbjUdpDoq7d3A2oJ2paxbY4XgYlYB2T/xsDB8HzgN9P/5V82YXzgM95WubyIGAHi8mcthvWY0Qt/8+Xp+pOQLt3sH5gHam7tTNpREojEtAe0dsdpXZx2vgnUiH+nVudvjjoXciLU9fFkNX7Dxdzu54mtLpT44/OHcH64GHriegMXXHSK/sA4WRCWjnnKHTJZePm9F731o37t6sOtHodHXm1nvhj8fE//ifb629AEfHG9fHzI+nLf35v4/9G7yDupz1Q3ffC7//7vRmem/lhHEJaPes9X2+Tq2pvj1sAzfr177ZetG6vP26eci7qtn3eqdAdxh4TOj76vhI1U++tZveu4Pjb5sP3dzH0Ji6/SPYAQpje8mAPq2VtwzBnAgoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIRE6o2EAAAZdSURBVAEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCAkoQEhAAUICChASUICQgAKEBBQgJKAAIQEFCP1/AJ5W750ZOacAAAAASUVORK5CYII=)