# Confidence Interval with Wilcoxon Test in Python for log-normal Distribution

I have a couple arrays of dN/dS scores, and I would like to calculate the confidence interval for each array of data. dN/dS scores are not normally distributed but are log-normally distrbuted, so I want some non-parametric approach for finding confidence intervals.

An answer in the CrossValidated forum suggested that in R you can use a Wilcoxon test to find the median-based confidence interval, a function that I do not think the python scipy.stats.wilcoxon() package does. Is there anyway of doing something similar in Python for the confidence intervals of log-normally distributed data?

For a example short array of dN/dS scores, here is the dN/dS scores of ATPases:

array([0.02288081, 0.44170839, 0.10549733, 0.17515196, 0.09449279,
0.07110412, 0.00893079, 0.23485109, 0.14533192, 0.05449631,
0.10281173, 0.02113355, 0.05157087, 0.01705113, 0.02651948,
0.09352947, 0.05828018, 0.20528157, 0.02873843, 0.0141598 ,
0.11262881, 0.06337332, 0.13689815, 0.10557703, 0.04452507,
0.05046484, 0.40241456, 0.05939199, 0.24423569, 0.05042892,
0.164257  , 0.03408215, 0.04737262, 0.01037463, 0.01246448,
0.02170375, 0.0241773 , 0.05136936, 0.02393366, 0.18913979,
0.15781334, 0.06448557, 0.04355384, 0.02821125, 0.08015629,
0.10985432, 0.06074574, 0.15775976, 0.05678278, 0.03749782,
0.05518756, 0.00770479, 0.21248167, 0.08005044, 0.16307954,
0.05783565, 0.05907416, 0.07044622, 0.13227131, 0.01627556,
0.10859962, 0.08149819, 0.05600647, 0.16098728, 0.15183062,
0.05202344, 0.01769589, 0.00789287, 0.07777749, 0.1324942 ,
0.12734709, 0.17146938, 0.03890857, 0.1296019 , 0.085146  ,
0.05602965, 0.02708089, 0.11807038, 0.07848828, 0.03291032,
0.36302533, 0.07800343, 0.06551307, 0.04676282, 0.04765273,
0.08060882, 0.06339636, 0.03349833, 0.01224308, 0.1481316 ,
0.31738452, 0.15690855, 0.0693822 , 0.020425  , 0.02909208,
0.03499225, 0.03019904, 0.13722717, 0.14403507, 0.01257245,
0.02223452, 0.07068784, 0.0544813 , 0.08738558, 0.02884046,
0.10549474, 0.06695546, 0.01341142, 0.09440411, 0.11840834,
0.08558889, 0.2688645 , 0.28313546, 0.15127967, 0.01463191,
0.1728421 ])

• I am not at sure about that, I would be amazed if SciPy couldn't do this. However in R is a base function. How many columns, comparisons? What are you measuring?
– M__
Oct 13, 2019 at 12:15

I guess you mean confidence interval of the mean (or median). I looked at the array you provided, even after log it does not look normal-distributed. Use the bootstrap for non-parametric values. I modified the code byJames Triveri, see example below, and you can use the mean or median. The textbook by Efron is a good reference for the method, which is also used in R package boot. There's also awebsite which should be ok.


import numpy as np

def bootstrap(data, n=1000, func=np.mean,p=0.95):

sample_size = len(data)
simulations = [func(np.random.choice(data, size=sample_size, replace=True)) for i in range(n)]
simulations.sort()
u_pval = (1+p)/2.
l_pval = (1-u_pval)
l_indx = int(np.floor(n*l_pval))
u_indx = int(np.floor(n*u_pval))
return(simulations[l_indx],simulations[u_indx])

x = np.array([0.02288081, 0.44170839, 0.10549733, 0.17515196, 0.09449279,
0.07110412, 0.00893079, 0.23485109, 0.14533192, 0.05449631,
0.10281173, 0.02113355, 0.05157087, 0.01705113, 0.02651948,
0.09352947, 0.05828018, 0.20528157, 0.02873843, 0.0141598 ,
0.11262881, 0.06337332, 0.13689815, 0.10557703, 0.04452507,
0.05046484, 0.40241456, 0.05939199, 0.24423569, 0.05042892,
0.164257  , 0.03408215, 0.04737262, 0.01037463, 0.01246448,
0.02170375, 0.0241773 , 0.05136936, 0.02393366, 0.18913979,
0.15781334, 0.06448557, 0.04355384, 0.02821125, 0.08015629,
0.10985432, 0.06074574, 0.15775976, 0.05678278, 0.03749782,
0.05518756, 0.00770479, 0.21248167, 0.08005044, 0.16307954,
0.05783565, 0.05907416, 0.07044622, 0.13227131, 0.01627556,
0.10859962, 0.08149819, 0.05600647, 0.16098728, 0.15183062,
0.05202344, 0.01769589, 0.00789287, 0.07777749, 0.1324942 ,
0.12734709, 0.17146938, 0.03890857, 0.1296019 , 0.085146  ,
0.05602965, 0.02708089, 0.11807038, 0.07848828, 0.03291032,
0.36302533, 0.07800343, 0.06551307, 0.04676282, 0.04765273,
0.08060882, 0.06339636, 0.03349833, 0.01224308, 0.1481316 ,
0.31738452, 0.15690855, 0.0693822 , 0.020425  , 0.02909208,
0.03499225, 0.03019904, 0.13722717, 0.14403507, 0.01257245,
0.02223452, 0.07068784, 0.0544813 , 0.08738558, 0.02884046,
0.10549474, 0.06695546, 0.01341142, 0.09440411, 0.11840834,
0.08558889, 0.2688645 , 0.28313546, 0.15127967, 0.01463191,
0.1728421 ])

# for mean,95% confidence interval, do
bootstrap(x,1000,np.mean,0.95)
# for median, 95% confidence interval, do
bootstrap(x,1000,np.median,0.95)


• This is a good answer and it would certainly work. This is the 95% CI of the mean and it's not the Wilcoxan. If the OP is using W it suggests the data is separated into columns and they are trying to determine the difference because different groups.
– M__
Oct 13, 2019 at 16:14
• Essentially the question isn't clear enough at present
– M__
Oct 13, 2019 at 16:15
• I agree with you @MichaelG. You can replace np.mean with np.median, its the same concept. Yes, let's hope VeritatemAmo provides more info. Oct 13, 2019 at 16:42