# Kubke Lab:Research/ABR/Notebook/2013/10/29

(Difference between revisions)
 Revision as of 20:40, 28 October 2013 (view source)m ← Previous diff Current revision (22:33, 28 October 2013) (view source) (→Feedback from Andy) (2 intermediate revisions not shown.) Line 11: Line 11: =General Entries= =General Entries= * Played around with filter() and butter() in library.signal [[User:M Fabiana Kubke|MF Kubke]] 20:11, 28 October 2013 (EDT) * Played around with filter() and butter() in library.signal [[User:M Fabiana Kubke|MF Kubke]] 20:11, 28 October 2013 (EDT) + + == Feedback from Andy == + + * Filter: Higher order is sharper cutoff. A sharp butterworth would be a 5th or 6th order.  I wouldn't use any higher. + + * Try an fft to see where the low pass cutoff should be. If the filter order is too high you get ringing. [ see http://t.co/QfrGaDjvj2 ] + + * entered wrong sampling rate - sampling rate is actually 50 kHz - need to recalculate the filters based on that. They are off by a factor of 10, so what says lowpassed at 100 is actually at 1k. =Personal Entries= =Personal Entries=

## Current revision

Hearing development in barn owls Main project page
Previous entry      Next entry

# General Entries

• Played around with filter() and butter() in library.signal MF Kubke 20:11, 28 October 2013 (EDT)

## Feedback from Andy

• Filter: Higher order is sharper cutoff. A sharp butterworth would be a 5th or 6th order. I wouldn't use any higher.
• Try an fft to see where the low pass cutoff should be. If the filter order is too high you get ringing. [ see http://t.co/QfrGaDjvj2 ]
• entered wrong sampling rate - sampling rate is actually 50 kHz - need to recalculate the filters based on that. They are off by a factor of 10, so what says lowpassed at 100 is actually at 1k.

# Personal Entries

## Fabiana

• File 2013-10-29-MFK.Rmd

# Trying to apply a low pass filter to ABR data

Using package signal

``````install.packages("signal")
``````
``````## Error: trying to use CRAN without setting a mirror
``````
``````library(signal)
``````
``````## Warning: package 'signal' was built under R version 3.0.2
``````
``````## Loading required package: MASS
##
## Attaching package: 'signal'
##
## The following object is masked from 'package:stats':
##
##     filter, poly
``````

``````dir()
``````
``````##  [1] "2013-10-27-MFK.Rmd"  "2013-10-27.html"     "2013-10-27.R"
##  [4] "2013-10-27.txt"      "2013-10-28-MFK.html" "2013-10-28-MFK.md"
##  [7] "2013-10-28-MFK.Rmd"  "2013-10-29-MFK.Rmd"  "Copy189L0A.ABR"
## [10] "Copy224l0a.txt"      "Copy233L0B.ABR.log"  "Copy233L0B.ABR.txt"
## [13] "figure"              "knitr-testr.html"    "knitr-testr.md"
## [16] "knitr-testr.Rmd"     "kntR2.html"          "kntR2.md"
## [19] "kntR2.Rmd"           "my first knitr.html" "my first knitr.md"
## [22] "my first knitr.Rmd"  "mydata_new"          "mydata_new.txt"
## [28] "signal.pdf"          "texput.log"
``````
``````mydata <- read.table("Copy233L0B.ABR.txt", header = T)
``````
``````##   msec left right spont  bin
## 1 0.00 1000  1000  1000 1000
## 2 0.02 1000  1000  1000 1000
## 3 0.04 1000  1000  1000 1000
## 4 0.06 1000  1000  1000 1000
## 5 0.08 1000  1000  1000 1000
## 6 0.10 1000  1000  1000 1000
``````
``````write.table(mydata, "mydata_new.txt")
``````
``````##   msec left right spont  bin
## 1 0.00 1000  1000  1000 1000
## 2 0.02 1000  1000  1000 1000
## 3 0.04 1000  1000  1000 1000
## 4 0.06 1000  1000  1000 1000
## 5 0.08 1000  1000  1000 1000
## 6 0.10 1000  1000  1000 1000
``````
``````mydata_new\$bincomp[1:1000] <- (mydata_new\$bin[1:1000] - (mydata_new\$left[1:1000] +
mydata_new\$right[1:1000]))
``````
``````##   msec left right spont  bin bincomp
## 1 0.00 1000  1000  1000 1000   -1000
## 2 0.02 1000  1000  1000 1000   -1000
## 3 0.04 1000  1000  1000 1000   -1000
## 4 0.06 1000  1000  1000 1000   -1000
## 5 0.08 1000  1000  1000 1000   -1000
## 6 0.10 1000  1000  1000 1000   -1000
``````
``````write.table(mydata_new, "mydata_new.txt")  #passed test
``````
``````##   msec left right spont  bin bincomp
## 1 0.00 1000  1000  1000 1000   -1000
## 2 0.02 1000  1000  1000 1000   -1000
## 3 0.04 1000  1000  1000 1000   -1000
## 4 0.06 1000  1000  1000 1000   -1000
## 5 0.08 1000  1000  1000 1000   -1000
## 6 0.10 1000  1000  1000 1000   -1000
``````

Test plotting for comparison

``````
plot(mydata_new[c(1, 2)], type = "l", main = "original left")
``````

Now - try to apply a low pass filterusing butter() to determine the properties of the filter.

``````bf <- butter(3, 0.1)
test1 <- filter(bf, mydata_new\$left)

par(mfrow = c(1, 2))
plot(mydata_new[c(1, 2)], type = "l", main = "original left")
plot(test1, main = "butter(3, 0.1)")
``````

Test type = “low” vs type = “high”

``````bf <- butter(3, 0.1, type = "low")
test1 <- filter(bf, mydata_new\$left)

par(mfrow = c(1, 2))
plot(test1, main = "butter type low")

bf <- butter(3, 0.1, type = "high")
test1 <- filter(bf, mydata_new\$left)
plot(test1, main = "butter type high")
``````

(guess default is low then)

Try changing filter properties - change filter order (first argument in function)

``````par(mfrow = c(2, 2))
plot(mydata_new[c(1, 2)], type = "l", main = "original left")

bf <- butter(1, 0.1, type = "low")
test1 <- filter(bf, mydata_new\$left)
plot(test1, main = "butter(1, 0.1)")

bf <- butter(3, 0.1, type = "low")
test1 <- filter(bf, mydata_new\$left)
plot(test1, main = "butter(3, 0.1)")

bf <- butter(15, 0.1, type = "low")
test1 <- filter(bf, mydata_new\$left)
plot(test1, main = "butter(30, 0.1)")
``````

Not really sure what the order is doing - ????

Stick with filter order == 3, change crirical frequency W. W-> [0-1] where 1 is the Nyquist frequency (½ of the sampling rate )

``````head(mydata_new)
``````
``````##   msec left right spont  bin bincomp
## 1 0.00 1000  1000  1000 1000   -1000
## 2 0.02 1000  1000  1000 1000   -1000
## 3 0.04 1000  1000  1000 1000   -1000
## 4 0.06 1000  1000  1000 1000   -1000
## 5 0.08 1000  1000  1000 1000   -1000
## 6 0.10 1000  1000  1000 1000   -1000
``````

sampling rate => every 0.2 msec ()

Looks like sampling rate at 5 kHz => Nyquist = 2.5 KHz

``````par(mfrow = c(2, 2))
plot(mydata_new[c(1, 2)], type = "l", main = "original left")

bf <- butter(3, 0.01, type = "low")
test1 <- filter(bf, mydata_new\$left)
plot(test1, main = "butter(3, 0.01 (25 Hz)")

bf <- butter(3, 0.1, type = "low")
test1 <- filter(bf, mydata_new\$left)
plot(test1, main = "butter(3, 0.1 (250 Hz))")

bf <- butter(3, 1, type = "low")
test1 <- filter(bf, mydata_new\$left)
plot(test1, main = "butter(30, 1 (2.5KHz))")
``````

Try the first derivative on the 250 LP

``````
bf <- butter(3, 0.1, type = "low")
mydata_new\$leftlow <- filter(bf, mydata_new\$left)

par(mfrow = c(1, 3))
plot(mydata_new\$msec, mydata_new\$left, type = "l", main = "original left")
plot(mydata_new\$msec, mydata_new\$leftlow, type = "l", main = "lowpassed left")
mydata_new\$diffleft2[2:1000] <- diff(mydata_new\$leftlow)/diff(mydata_new\$msec)
plot(mydata_new\$msec, mydata_new\$diffleft2, ylim = c(-3000, 3000), type = "l",
main = "first diff on lowpassed left")
``````

First deriv stil looks noisy - try different filters

``````
bf <- butter(3, 0.04, type = "low")
mydata_new\$leftlow <- filter(bf, mydata_new\$left)

par(mfrow = c(1, 3))
plot(mydata_new\$msec, mydata_new\$left, type = "l", main = "original left")
plot(mydata_new\$msec, mydata_new\$leftlow, type = "l", main = "lowpassed  100 left")
mydata_new\$diffleft2[2:1000] <- diff(mydata_new\$leftlow)/diff(mydata_new\$msec)
plot(mydata_new\$msec, mydata_new\$diffleft2, ylim = c(-3000, 3000), type = "l",
main = "first diff on lowpassed left")
``````

``````
par(mfrow = c(1, 1))

plot(mydata_new\$msec, mydata_new\$diffleft2, ylim = c(-3000, 3000), type = "l",
main = "first diff on lowpassed left")
``````

Mmmmmmm….

## Andy

• Enter content here

## Oris

• Enter content here