Sparse PLS on R

# Reference

* [An Introduction to the ‘spls’ Package, Version 1.0](https://cran.r-project.org/web/packages/spls/vignettes/spls-example.pdf)
* [Sparse Partial Least Squares Regression for Simultaneous Dimension Reduction and Variable Selection](http://pages.stat.wisc.edu/~keles/Papers/SPLS_Nov07.pdf)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/278468/8e213f85-496b-accc-47fc-c672ece7e02a.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/278468/ff2240ed-d1ef-36c0-2943-a9a46641a732.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/278468/fce11d69-fbf3-b1ee-29c0-09869279e76b.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/278468/2ddc3a6e-d76c-c9d1-4210-765c275748c1.png)

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/278468/f0fd38de-d8be-1483-315b-6959f2a3d4ef.png)

# library


```R
install.packages("spls")
```

    Updating HTML index of packages in '.Library'
    
    Making 'packages.html' ...
     done
    



```R
library(spls)
```

    Sparse Partial Least Squares (SPLS) Regression and
    Classification (version 2.2-3)
    
    



```R
library(ggplot2)

```

# data


```R
data(yeast)
```

xはクロマチンの免疫沈降、yは発現サイクル？


```R
help(yeast)
```



<table width="100%" summary="page for yeast {spls}"><tr><td>yeast {spls}</td><td style="text-align: right;">R Documentation</td></tr></table>

<h2>Yeast Cell Cycle Dataset</h2>

<h3>Description</h3>

<p>This is the Yeast Cell Cycle dataset used in Chun and Keles (2010).
</p>


<h3>Usage</h3>

<pre> data(yeast) </pre>


<h3>Format</h3>

<p>A list with two components:
</p>

<dl>
<dt>x</dt><dd><p> ChIP-chip data. A matrix with 542 rows and 106 columns.</p>
</dd>
<dt>y</dt><dd><p> Cell cycle gene expression data.
A matrix with 542 rows and 18 columns.</p>
</dd>
</dl>



<h3>Details</h3>

<p>Matrix <code>y</code> is cell cycle gene expression data (Spellman et al., 1998)
of 542 genes from an <i>&alpha;</i> factor based experiment.
Each column corresponds to mRNA levels
measured at every 7 minutes during 119 minutes (a total of 18 measurements).
Matrix <code>x</code> is the chromatin immunoprecipitation on chip (ChIP-chip) data of
Lee et al. (2002) and it contains the binding information for 106
transcription factors. See Chun and Keles (2010) for more details. 
</p>


<h3>Source</h3>

<p>Lee TI, Rinaldi NJ, Robert F, Odom DT, Bar-Joseph Z, Gerber GK, Hannett NM,
Harbison CT, Thomson CM, Simon I, Zeitlinger J, Jennings EG, Murray HL,
Gordon DB, Ren B, Wyrick JJ, Tagne JB, Volkert TL, Fraenkel E, Gifford DK,
and Young RA (2002), &quot;Transcriptional regulatory networks in <em>Saccharomyces cerevisiae</em>&quot;,
<em>Science</em>, Vol. 298, pp. 799&ndash;804.
</p>
<p>Spellman PT, Sherlock G, Zhang MQ, Iyer VR, Anders K, Eisen MB, Brown PO,
Botstein D, and Futcher B (1998), &quot;Comprehensive identification of cell cycle-regulated genes of
the yeast <em>Saccharomyces cerevisiae</em> by microarray hydrization&quot;,
<em>Molecular Biology of the Cell</em>, Vol. 9, pp. 3273&ndash;3279.
</p>


<h3>References</h3>

<p>Chun H and Keles S (2010), &quot;Sparse partial least squares
for simultaneous dimension reduction and variable selection&quot;,
<em>Journal of the Royal Statistical Society - Series B</em>, Vol. 72, pp. 3&ndash;25. </p>


<h3>Examples</h3>

<pre>
data(yeast)
yeast$x[1:5,1:5]
yeast$y[1:5,1:5]
</pre>

<hr /><div style="text-align: center;">[Package <em>spls</em> version 2.2-3 ]</div>



```R
yeast$x[1:5,1:5]
```


<table>
<caption>A matrix: 5 × 5 of type dbl</caption>
<thead>
	<tr><th></th><th scope=col>ABF1_YPD</th><th scope=col>ACE2_YPD</th><th scope=col>ADR1_YPD</th><th scope=col>ARG80_YPD</th><th scope=col>ARG81_YPD</th></tr>
</thead>
<tbody>
	<tr><th scope=row>21</th><td>-0.2722730</td><td> 0.21932294</td><td> 0.9238359567</td><td>-0.4755756</td><td>-0.10389318</td></tr>
	<tr><th scope=row>41</th><td> 0.1691280</td><td> 0.53831198</td><td> 0.0097604993</td><td>-0.3219534</td><td>-0.19750606</td></tr>
	<tr><th scope=row>71</th><td>-0.1388962</td><td> 0.02636382</td><td> 0.0877516229</td><td>-0.2234093</td><td> 0.10307741</td></tr>
	<tr><th scope=row>78</th><td>-0.2865169</td><td>-0.31409427</td><td>-0.0454998435</td><td> 0.3262217</td><td> 0.27757502</td></tr>
	<tr><th scope=row>102</th><td>-0.4950561</td><td>-0.14827419</td><td> 0.0002987512</td><td>-0.2179458</td><td>-0.02539585</td></tr>
</tbody>
</table>




```R
plot(yeast$x[,1])
```


![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/278468/a16cc5b4-0779-4043-513d-09ae12e8d5bf.png)




```R
yeast$y[1:5,1:5]
```


<table>
<caption>A matrix: 5 × 5 of type dbl</caption>
<thead>
	<tr><th></th><th scope=col>alpha0</th><th scope=col>alpha7</th><th scope=col>alpha14</th><th scope=col>alpha21</th><th scope=col>alpha28</th></tr>
</thead>
<tbody>
	<tr><th scope=row>1</th><td>-0.36</td><td>-0.42</td><td> 0.29</td><td>-0.14</td><td>-0.19</td></tr>
	<tr><th scope=row>2</th><td> 1.04</td><td> 0.19</td><td> 0.47</td><td>-1.03</td><td>-0.63</td></tr>
	<tr><th scope=row>5</th><td>-0.30</td><td>-0.45</td><td> 0.75</td><td> 0.37</td><td> 0.27</td></tr>
	<tr><th scope=row>8</th><td>-0.46</td><td> 0.12</td><td>-0.06</td><td>-0.76</td><td>-0.70</td></tr>
	<tr><th scope=row>9</th><td>-1.35</td><td>-0.86</td><td>-0.22</td><td>-0.38</td><td>-0.65</td></tr>
</tbody>
</table>




```R
plot(yeast$y[,1])
```


![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/278468/8806cd96-0770-5af2-f61a-0a1c627af68a.png)



↓ 変数が542個もある・・


```R
str(yeast)
```

    List of 2
     $ x: num [1:542, 1:106] -0.272 0.169 -0.139 -0.287 -0.495 ...
      ..- attr(*, "dimnames")=List of 2
      .. ..$ : chr [1:542] "21" "41" "71" "78" ...
      .. ..$ : chr [1:106] "ABF1_YPD" "ACE2_YPD" "ADR1_YPD" "ARG80_YPD" ...
     $ y: num [1:542, 1:18] -0.36 1.04 -0.3 -0.46 -1.35 -2.06 -1.61 -0.07 0.11 0.15 ...
      ..- attr(*, "dimnames")=List of 2
      .. ..$ : chr [1:542] "1" "2" "5" "8" ...
      .. ..$ : chr [1:18] "alpha0" "alpha7" "alpha14" "alpha21" ...

# SPLS

```R
set.seed(1)
cv <- cv.spls(yeast$x, yeast$y, eta = seq(0.1,0.9,0.1), K = c(5:10))
```

    eta = 0.1 
    eta = 0.2 
    eta = 0.3 
    eta = 0.4 
    eta = 0.5 
    eta = 0.6 
    eta = 0.7 
    eta = 0.8 
    eta = 0.9 
    
    Optimal parameters: eta = 0.6, K = 8



![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/278468/caa827d7-0346-af6d-7565-32e735743f27.png)




```R
f <- spls(yeast$x, yeast$y, eta = cv$eta.opt, K = cv$K.opt)
```


```R
print(f)
```

    
    Sparse Partial Least Squares for multivariate responses
    ----
    Parameters: eta = 0.6, K = 8, kappa = 0.5
    PLS algorithm:
    pls2 for variable selection, simpls for model fitting
    
    SPLS chose 56 variables among 106 variables
    
    Selected variables: 
    ACE2_YPD	ARG80_YPD	ARG81_YPD	ASH1_YPD	AZF1_YPD	
    BAS1_YPD	CBF1_YPD	CHA4_YPD	CRZ1_YPD	FHL1_YPD	
    FKH1_YPD	FKH2_YPD	FZF1_YPD	GAT1_YPD	GAT3_YPD	
    GCN4_YPD	GCR2_YPD	GLN3_YPD	HAA1_YPD	HAP2_YPD	
    HAP5_YPD	HIR1_YPD	HIR2_YPD	IME4_YPD	INO4_YPD	
    A1..MATA1._YPD	MBP1_YPD	MCM1_YPD	MET4_YPD	MSN2_YPD	
    NDD1_YPD	NRG1_YPD	PHD1_YPD	PHO2_YPD	PUT3_YPD	
    RCS1_YPD	REB1_YPD	RFX1_YPD	RIM101_YPD	RME1_YPD	
    RTG1_YPD	RTG3_YPD	SIP4_YPD	SOK2_YPD	STB1_YPD	
    STE12_YPD	STP2_YPD	SWI4_YPD	SWI5_YPD	SWI6_YPD	
    THI2_YPD	YAP1_YPD	YAP6_YPD	YAP7_YPD	YFL044C_YPD	
    YJL206C_YPD	



```R
coef.f <- coef(f)
```


```R
coef.f[1:5,1:5]
```


<table>
<caption>A matrix: 5 × 5 of type dbl</caption>
<thead>
	<tr><th></th><th scope=col>alpha0</th><th scope=col>alpha7</th><th scope=col>alpha14</th><th scope=col>alpha21</th><th scope=col>alpha28</th></tr>
</thead>
<tbody>
	<tr><th scope=row>ABF1_YPD</th><td> 0.0000000</td><td> 0.000000000</td><td>0.00000000</td><td> 0.0000000000</td><td> 0.000000000</td></tr>
	<tr><th scope=row>ACE2_YPD</th><td> 0.0874325</td><td> 0.068452293</td><td>0.01374781</td><td>-0.0002541969</td><td>-0.033302624</td></tr>
	<tr><th scope=row>ADR1_YPD</th><td> 0.0000000</td><td> 0.000000000</td><td>0.00000000</td><td> 0.0000000000</td><td> 0.000000000</td></tr>
	<tr><th scope=row>ARG80_YPD</th><td>-0.0486881</td><td>-0.019092797</td><td>0.02063442</td><td> 0.0300421634</td><td> 0.007925553</td></tr>
	<tr><th scope=row>ARG81_YPD</th><td>-0.0168849</td><td> 0.009465868</td><td>0.06353825</td><td> 0.0541704059</td><td> 0.006978985</td></tr>
</tbody>
</table>




```R
plot.spls(f, yvar=1 )
```


![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/278468/07f8f4ff-d4b0-1146-8d51-aed99ddd9829.png)




```R
coefplot.spls( f, nwin=c(2,2), xvar=c(1:4) )
```


![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/278468/8c9053a9-820b-80e6-1e03-9aeba8c11cbd.png)