## Tuesday, January 20, 2015

### Adjustment for Multiple Comparison Tests with R: Resources on the web

1. Bonferroni correction

p.adjust(p, method = "bonferroni")

2. Sidak (Dunn-Sidak) correction

3. Holm-Bonferroni correction

p.adjust(p, method = "holm")

4. Hochberg correction

p.adjust(p, method = "hochberg")

5. Hommel correction

p.adjust(p, method = "hommel")

6. Benjamini-Hochberg correction

p.adjust(p, method = "BH")
or equivalently
p.adjust(p, method = "fdr")

7. Benjamini–Yekutieli (Benjamini–Hochberg–Yekutieli) correction

p.adjust(p, method = "BY")

## Wednesday, February 5, 2014

### ggPlot2: Histogram with jittered stripchart

Here is an example of a Histogram plot, with a stripchart (vertically jittered) along the x side of the plot.

Alternatively, using the geom_rug function:

Of course this simplicistic method need to be adjusted in vertical position of the stripchart or rugchart (y=-2, here), and the relative proportion of points jittering. Continue reading »

## Sunday, February 2, 2014

### Boxplot with mean and standard deviation in ggPlot2 (plus Jitter)

When you create a boxplot in R, it automatically computes median, first and third quartile ("hinges") and 95% confidence interval of median ("notches").

But we would like to change the default values of boxplot graphics with the mean, the mean + standard deviation, the mean - S.D., the min and the max values.
Here is an example solved using ggplot2 package. Plus here are represented points (the single values) jittered horizontally.

## Saturday, September 17, 2011

### Implementation of the CDC Growth Charts in R

I implemented in R a function to re-create the CDC Growth Chart, according to the data provided by the CDC.

In order to use this function, you need to download the .rar file available at this megaupload link.

Then unrar the file, and put the Growth folder in your main directory, as selected in R. You are now able to use the two functions i'm going to illustrate.

growthFun.R

The function `growthFun` allows you to draw 8 different growth chart, which are different for Male and Female (sixteen in total).
The only parameters you need to input are:
` sex = c("m", "f")`
` type = c("wac36", "lac36", "wlc", "hac", "wsc", "wac20", "lac20", "bac")`
The explanation for the type's parameters code are in the first part of the function code.
Eventually you can modify the `pat` variable, if you want to put the Growth folder in another place (not in the main directory of R).

I reccomend to use the `pdf()` graphic device for best resolution.

Hese is an example of the output you can obtain, with the following code:

`pdf("hac_example.pdf", paper="a4", width=0, height=0)growthFun("m", "hac")dev.off()`

MygrowthFun.R

The function `MygrowthFun` allows you to personalize the output of the previous function, with specific patient's data.
The parameters you can modify are:
` sex=c("m", "f")type=c("wac36", "lac36", "wlc", "hac", "wsc", "wac20", "lac20", "bac", "bmi.adv")path="./Growth/"name = NULLsurname = NULLbirth_date = NULLmydataAA = NULL`

The three parameter `sex`, `type` and `path` are the same of the `growthFun` function. The three parameters `name`, `surname` and `birth_date` refer to the patient's data; you can add this data in form of `character()`.
`mydataAA` is an optional parameters with the values measured on your patients during the time you follow up him. Generally you need to input this data in form of a `data.frame()`.
In the `type` parameter there is an additional choice: `bmi.adv` allows you to obtain three chart (`wac20, lac20, bac` - see the explanation codes), if your `mydataAA` dataframe contains data about Stature and Weight during the time of follow up.

Details.
Let's see the format of `mydataAA`, according to the `type` of chart you want to graph.

`type = wac36mydataAA: first column = months of measurement, from 0 to 36second column = weight (in kg)type = lac36mydataAA: first column = months of measurement, from 0 to 36second column = length (in cm)type = hacmydataAA: first column = months of measurement, from 0 to 36second column = head circumference (in cm)type = wac20mydataAA: first column = months of measurement, from 24 to 240 (from 2 to 20 years)second column = weight (in kg)type = lac20mydataAA: first column = months of measurement, from 24 to 240 (from 2 to 20 years)second column = stature (in cm)type = bmi.advmydataAA: first column (months) = months of measurement, from 24 to 240 (from 2 to 20 years)second column (stature) = stature (in cm)third column (weight)= weight (in kg)`

In the last type it's not importat the order of the columns, but here are important their names.

Examples.
Let's see some example. Suppose that you are following the growth of a new born (her name is Alyssa Gigave, born on 07/08/2009), and you collect the following data:

`Months  Length0       502       553       565       618       719       7212      7515      7518      8121      8926      9127      9430      9535      98`

So you can create your personalized graph in this way:

`alyssa_data <- data.frame(   months=c(0, 2, 3, 5, 8, 9, 12, 15, 18, 21, 26, 27, 30, 35),   length=c(50, 55, 56, 61, 71, 72, 75, 75, 81, 89, 91, 94, 95, 98))  pdf("alyssa_growth_chart.pdf", paper="a4", width=0, height=0)  MygrowthFun(sex="f", type="lac36", name="Alyssa", surname="Gigave", birth_date="july 08, 2009", mydataAA=alyssa_data)  dev.off()`

The output is the following pdf file:

Now suppose that you're a pediatric doctor, and that you follow a boy (Tommy Cigalino, born on 07/08/1980). Whenever he has come to you, you collect his weight and stature, and the months from his birth he was. So you have the following data:

`  months stature weight     25      98     17     31     100     21     34     102     27     35     104     29     58     106     30     60     109     32     70     111     33     85     118     34     88     119     36     89     120     39     91     121     42    102     126     45    107     128     47    108     135     49    120     144     51    134     145     52    154     148     54    166     152     55    169     157     62    170     158     63    178     163     64    179     167     68    181     168     71    219     169     74    234     176     76`

So you can create three graphs (`wac20, lac20, bac`), using the `bmi.adv` type:

`tommy_data <- data.frame(  months = c( 25, 31, 34, 35, 58, 60,               70, 85, 88, 89, 91, 102,               107, 108, 120, 134, 154,               166, 169, 170, 178, 179,               181, 219, 234),   stature = c( 98, 100, 102, 104, 106,               109, 111, 118, 119, 120,               121, 126, 128, 135, 144,               145, 148, 152, 157, 158,               163, 167, 168, 169, 176),   weight = c( 17, 21, 27, 29, 30, 32,               33, 34, 36, 39, 42, 45,               47, 49, 51, 52, 54, 55,               62, 63, 64, 68, 71, 74,               76))  pdf("tommy_growth_chart.pdf", paper="a4", width=0, height=0)  MygrowthFun(sex="m", type="bmi.adv", name="Tommy", surname="Cigalino", birth_date="july 08, 1980", mydataAA=tommy_data)  dev.off()`

Tommaso MARTINO, 17/09/2011

REFERENCES

• http://www.cdc.gov/growthcharts/cdc_charts.htm

• http://www.cdc.gov/growthcharts/clinical_charts.htm

• http://www.cdc.gov/growthcharts/percentile_data_files.htm

• Kuczmarski RJ, Ogden CL, Guo, SS, et al. CDC growth charts for the United States: Methods and Development. Vital Health Stat; 11 (246) National Center for Health Statistics. 2002.

## Wednesday, September 7, 2011

### R is a cool sound editor!

Capabilities of R are definitely unless! After my previous posts about some easy image editing in R (they are here, and here), now is the time to explore if R is capable of sound editing!

Just for fun, here I created a function that receives a phone number (or another sequence of numbers), and returns the equivalent melody you can listen if you press that sequence on your house' phone... =D

It requires the `sound` library, and here's the code.

Now you can simply create your phone melody =)

`s2 <- PlayTel("556c885a4623#")`

You can listen to it with the command:

`play(s2)`

(NOTE: in Windows 7 I was unable to find a wave player that works on batch mode - i.e. mplay32.exe. So this command doesn't work on Windows 7. It works on Windows XP)

You can save the output using the command:

`saveSample(s2, "tel.wav")`

(This command works on Windows 7)

Here is an example of the output:

Have fun!! =)