• API
• FAQ
• Tools
• Archive
daily pastebin goal
49%
SHARE
TWEET

# Untitled

a guest Feb 23rd, 2018 64 Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
1. #Pre-class work
2. dat=data.frame(x=c(1,2,3,4,5,6),
3.                y=c(1,3,5,6,8,12))
4.
5. min.RSS <- function(data, par) {
6.   with(data, sum((par[1] + par[2] * x - y)^2))
7. }
8.
9. result <- optim(par = c(0, 1), min.RSS, data = dat)
10.
11. library(Matching)
12. data(lalonde)
13.
15.
16. lm(re78~educ, data=lalonde)
17.
18. 't
19. Optimization:
20.
21. \$par
22. [1] 920.8695 429.5908
23.
24. \$value
25. [1] 19262169077
26.
27. \$counts
29.      123       NA
30.
31. \$convergence
32. [1] 0
33.
34. \$message
35. NULL
36.
37. Regression:
38.
39. lm(formula = re78 ~ educ, data = lalonde)
40.
41. Coefficients:
42. (Intercept)         educ
43.       918.2        429.9
44.
45. 't
46.
47. claw <- function(xx) {
48.   x <- xx[1]
49.   y <- (0.46*(dnorm(x,-1.0,2.0/3.0) + dnorm(x,1.0,2.0/3.0)) +
50.           (1.0/300.0)*(dnorm(x,-0.5,.01) + dnorm(x,-1.0,.01) + dnorm(x,-1.5,.01)) +
51.           (7.0/300.0)*(dnorm(x,0.5,.07) + dnorm(x,1.0,.07) + dnorm(x,1.5,.07)))
52.   return(y)
53. }
54.
55. #Claw
56.
57. #clawx=sample((-200:200),replace=1)
58. #clawdata=data.frame(y=sapply(clawx,claw),x=clawx)
59. invclaw=function(x){-claw(x)}
60. optim(par = -2, invclaw)
61. optim(par = 0, invclaw)
62. optim(par = 2, invclaw)
63. #optim(par=2,claw,control=c(5,-1,0.01,1e-3,100,0,1e-8,0,10,1,5,1e7,0,10,10))
64. optimize(claw,c(-2,2),maximum=1)
65.
66. #Contaminated data
67. set.seed(123)
68. # Define the x values
69. x <- c(1:12)
70.
71. # Define the y values for the first 12 observations, in terms of x
72. # y = -x + 10 + error term
73. y <- c(-1*x[1:12] + 10 + rnorm(12, sd = 1))
74.
75. lm(y~x)
76.
77. min.reg <- function(par) {
78.   loss <- mean((y - (par[1]*x + par[2]))^2)
79.   return(loss)
80. }
81.
82. result.reg <- optim(par = c(1,1), min.reg)
83. result.reg\$par
84. result.reg\$value
85.
86. # Introduce a totally different data generating process for 1 observation
87. # (Maybe you have some small 'data contamination'): (21, 21)
88. x <- c(x,21)
89. y <- c(y, 21)
90.
91. 't
92. Use R to obtain the coefficients for the simple regression
93. Use optim to reproduce those coefficents
94. Plot the regression line
95. lm1 <- reg(y~x)
96. plot(x,y)
97. abline(lm1, col = “blue”, lwd = 3)
98. Consider what’s wrong with the regression line
99. Use optim() to run a robust regression (robust to data contamination) that minimizes the median squared residual instead of the mean of the squared residuals. Try it with the default optim() method, and then experiment with different starting values. Then try it with the “SANN” method, repeating the process with different starting values.
100. Using your best results, add a robust regression line to your plot using a different color than the one you used in (3) above.
101. t'
RAW Paste Data
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy.

Top