You signed in with another tab or window. Reload to refresh your session.You signed out in another tab or window. Reload to refresh your session.You switched accounts on another tab or window. Reload to refresh your session.Dismiss alert
Copy file name to clipboardExpand all lines: episodes/smoothing.Rmd
+52-34Lines changed: 52 additions & 34 deletions
Original file line number
Diff line number
Diff line change
@@ -5,12 +5,12 @@ exercises: 30
5
5
---
6
6
7
7
:::::: questions
8
-
-"Can a model be fitted to a dataset which shape is unknown but smooth?"
8
+
- Can a model be fitted to a dataset which shape is unknown but smooth?
9
9
::::::
10
10
11
11
:::::: objectives
12
-
-"Fit a smooth regression model to data which behavior depends conditionally on a set of predictors"
13
-
-"Predict the expected value of a smooth model given the value of the predictors"
12
+
- Fit a smooth regression model to data which behavior depends conditionally on a set of predictors.
13
+
- Predict the expected value of a smooth model given the value of the predictors.
14
14
::::::
15
15
16
16
@@ -19,20 +19,26 @@ exercises: 30
19
19
20
20
## Smoothing
21
21
22
-
Smoothing is a very powerful technique used all across data analysis. It is designed to estimate $f(x)$ when the shape is unknown, but assumed to be _smooth_. The general idea is to group data points that are expected to have similar expectations and compute the average, or fit a simple parametric model. We illustrate two smoothing techniques using a gene expression example.
22
+
Smoothing is a very powerful technique used all across data analysis. It is
23
+
designed to estimate $f(x)$ when the shape is unknown, but assumed to be
24
+
_smooth_. The general idea is to group data points that are expected to have
25
+
similar expectations and compute the average, or fit a simple parametric model.
26
+
We illustrate two smoothing techniques using a gene expression example.
23
27
24
28
The following data are gene expression measurements from replicated RNA samples.
25
29
26
30
27
31
```{r, message=FALSE, warning=FALSE}
28
-
##Following three packages are available from Bioconductor
32
+
##Following three packages are available from Bioconductor
29
33
library(Biobase)
30
34
library(SpikeIn)
31
35
library(hgu95acdf)
32
36
data(SpikeIn95)
33
37
```
34
38
35
-
We consider the data used in an MA-plot comparing two replicated samples ($Y$ = log ratios and $X$ = averages) and take down-sample in a way that balances the number of points for different strata of $X$ (code not shown):
39
+
We consider the data used in an MA-plot comparing two replicated samples
40
+
($Y$ = log ratios and $X$ = averages) and take down-sample in a way that
41
+
balances the number of points for different strata of $X$ (code not shown):
36
42
37
43
```{r, echo=FALSE}
38
44
## Example with two columns
@@ -64,9 +70,10 @@ mypar()
64
70
plot(X,Y)
65
71
```
66
72
67
-

68
-
69
-
In the MA plot we see that $Y$ depends on $X$. This dependence must be a bias because these are based on replicates, which means $Y$ should be 0 on average regardless of $X$. We want to predict $f(x)=\mbox{E}(Y \mid X=x)$ so that we can remove this bias. Linear regression does not capture the apparent curvature in $f(x)$:
73
+
In the MA plot we see that $Y$ depends on $X$. This dependence must be a bias
74
+
because these are based on replicates, which means $Y$ should be 0 on average regardless of $X$. We want to predict $f(x)=\mbox{E}(Y \mid X=x)$ so that we can
75
+
remove this bias. Linear regression does not capture the apparent curvature in
76
+
$f(x)$:
70
77
71
78
```{r MAplot_with_regression_line, fig.cap="MA-plot comparing gene expression from two arrays with fitted regression line. The two colors represent positive and negative residuals.",fig.width=10.5,fig.height=5.25}

80
-
81
86
The points above the fitted line (green) and those below (purple) are not evenly distributed. We therefore need an alternative more flexible approach.
82
87
83
88
## Bin Smoothing
84
89
85
-
Instead of fitting a line, let's go back to the idea of stratifying and computing the mean. This is referred to as _bin smoothing_. The general idea is that the underlying curve is "smooth" enough so that, in small bins, the curve is approximately constant. If we assume the curve is constant, then all the $Y$ in that bin have the same expected value. For example, in the plot below, we highlight points in a bin centered at 8.6, as well as the points of a bin centered at 12.1, if we use bins of size 1. We also show the fitted mean values for the $Y$ in those bins with dashed lines (code not shown):
90
+
Instead of fitting a line, let's go back to the idea of stratifying and
91
+
computing the mean. This is referred to as _bin smoothing_. The general idea is
92
+
that the underlying curve is "smooth" enough so that, in small bins, the curve
93
+
is approximately constant. If we assume the curve is constant, then all the $Y$
94
+
in that bin have the same expected value. For example, in the plot below, we
95
+
highlight points in a bin centered at 8.6, as well as the points of a bin
96
+
centered at 12.1, if we use bins of size 1. We also show the fitted mean values
97
+
for the $Y$ in those bins with dashed lines (code not shown):
86
98
87
99
```{r binsmoother, fig.cap="MAplot comparing gene expression from two arrays with bin smoother fit shown for two points.",fig.width=10.5,fig.height=5.25, echo=FALSE}

107
-
108
-
By computing this mean for bins around every point, we form an estimate of the underlying curve $f(x)$. Below we show the procedure happening as we move from the smallest value of $x$ to the largest. We show 10 intermediate cases as well (code not shown):
118
+
By computing this mean for bins around every point, we form an estimate of the underlying curve $f(x)$. Below we show the procedure happening as we move from
119
+
the smallest value of $x$ to the largest. We show 10 intermediate cases as well
120
+
(code not shown):
109
121
110
122
```{r bin_smoothing_demo, fig.cap="Illustration of how bin smoothing estimates a curve. Showing 12 steps of process.",fig.width=10.25,fig.height=10.25, echo=FALSE}
111
123
windowSize <- 0.5
@@ -125,8 +137,6 @@ for(i in seq(along=centers)){
125
137
}
126
138
```
127
139
128
-

129
-
130
140
The final result looks like this (code not shown):
131
141
132
142
```{r bin_smooth_final, fig.cap="MA-plot with curve obtained with bin-smoothed curve shown.", fig.width=10.5,fig.height=5.25,echo=FALSE}
@@ -135,14 +145,19 @@ plot(X, Y, col="darkgrey", pch=16)
135
145
lines(centers, smooth, col="black", lwd=3)
136
146
```
137
147
138
-

139
-
140
-
There are several functions in R that implement bin smoothers. One example is `ksmooth`. However, in practice, we typically prefer methods that use slightly more complicated models than fitting a constant. The final result above, for example, is somewhat wiggly. Methods such as `loess`, which we explain next, improve on this.
148
+
There are several functions in R that implement bin smoothers. One example is `ksmooth`. However, in practice, we typically prefer methods that use slightly
149
+
more complicated models than fitting a constant. The final result above, for
150
+
example, is somewhat wiggly. Methods such as `loess`, which we explain next,
151
+
improve on this.
141
152
142
153
## Loess
143
154
144
-
Local weighted regression (loess) is similar to bin smoothing in principle. The main difference is that we approximate the local behavior with a line or a parabola. This permits us to expand the bin sizes, which stabilizes the estimates. Below we see lines fitted to two bins that are slightly larger than those we used for the bin smoother (code not shown). We can use larger bins because fitting lines provide slightly more flexibility.
145
-
155
+
Local weighted regression (loess) is similar to bin smoothing in principle. The
156
+
main difference is that we approximate the local behavior with a line or a
157
+
parabola. This permits us to expand the bin sizes, which stabilizes the
158
+
estimates. Below we see lines fitted to two bins that are slightly larger than
159
+
those we used for the bin smoother (code not shown). We can use larger bins
160
+
because fitting lines provide slightly more flexibility.
146
161
147
162
```{r loess, fig.cap="MA-plot comparing gene expression from two arrays with bin local regression fit shown for two points.",fig.width=10.5,fig.height=5.25, echo=FALSE}
148
163
centers <- seq(min(X), max(X),0.1)
@@ -167,9 +182,8 @@ a <- min(X[ind]); b <- max(X[ind])

171
-
172
-
As we did for the bin smoother, we show 12 steps of the process that leads to a loess fit (code not shown):
185
+
As we did for the bin smoother, we show 12 steps of the process that leads to a
186
+
loess fit (code not shown):
173
187
174
188
```{r loess_demo, fig.cap="Illustration of how loess estimates a curve. Showing 12 steps of the process.",fig.width=10.25,fig.height=10.25, echo=FALSE}
175
189
mypar (4,3)
@@ -193,18 +207,15 @@ for(i in seq(along=centers)){
193
207
}
194
208
```
195
209
196
-

197
-
198
-
The final result is a smoother fit than the bin smoother since we use larger sample sizes to estimate our local parameters (code not shown):
210
+
The final result is a smoother fit than the bin smoother since we use larger
211
+
sample sizes to estimate our local parameters (code not shown):
199
212
200
213
```{r loess_final, fig.cap="MA-plot with curve obtained with loess.", fig.width=10.5,fig.height=5.25,echo=FALSE}
201
214
mypar (1,1)
202
215
plot(X,Y,col="darkgrey",pch=16)
203
216
lines(centers,smooth,col="black",lwd=3)
204
217
```
205
218
206
-

207
-
208
219
The function `loess` performs this analysis for us:
209
220
210
221
```{r loess2, fig.cap="Loess fitted with the loess function.", fig.width=10.25,fig.height=5.25}

222
-
223
-
There are three other important differences between `loess` and the typical bin smoother. The first is that rather than keeping the bin size the same, `loess` keeps the number of points used in the local fit the same. This number is controlled via the `span` argument which expects a proportion. For example, if `N` is the number of data points and `span=0.5`, then for a given $x$ , `loess` will use the `0.5*N` closest points to $x$ for the fit. The second difference is that, when fitting the parametric model to obtain $f(x)$, `loess` uses weighted least squares, with higher weights for points that are closer to $x$. The third difference is that `loess` has the option of fitting the local model robustly. An iterative algorithm is implemented in which, after fitting a model in one iteration, outliers are detected and downweighted for the next iteration. To use this option, we use the argument `family="symmetric"`.
232
+
There are three other important differences between `loess` and the typical bin smoother. The first is that rather than keeping the bin size the same, `loess`
233
+
keeps the number of points used in the local fit the same. This number is
234
+
controlled via the `span` argument which expects a proportion. For example, if
235
+
`N` is the number of data points and `span=0.5`, then for a given $x$ , `loess`
236
+
will use the `0.5*N` closest points to $x$ for the fit. The second difference is
237
+
that, when fitting the parametric model to obtain $f(x)$, `loess` uses weighted
238
+
least squares, with higher weights for points that are closer to $x$. The third difference is that `loess` has the option of fitting the local model robustly.
239
+
An iterative algorithm is implemented in which, after fitting a model in one
240
+
iteration, outliers are detected and downweighted for the next iteration. To use
241
+
this option, we use the argument `family="symmetric"`.
224
242
225
243
::::::::::::::::::::::::::::::::::::: challenge
226
244
@@ -323,5 +341,5 @@ popsd(res)
323
341
:::::::::::::::::::::::::::::::::::::::::::::::
324
342
325
343
:::::: keypoints
326
-
-"The smoothing methods work well when used inside the range of predictor values seen in the training set, however them are not suitable for extrapolation the prediction outside those ranges."
344
+
- The smoothing methods work well when used inside the range of predictor values seen in the training set, however them are not suitable for extrapolation the prediction outside those ranges.
0 commit comments