How to test bimodal data of two factors？ I have an unbalanced data frame. I want to check if the difference between treatment and fraction and their interaction can cause the “conne_density_pixel” to be different (especially treatment effect, because I can intuitively feel from the excel table that treatment effect may be more significant than fraction effect ).
This is my data for your reference.
#read data
input <- read.csv("input.csv",sep=",",header=TRUE)

#this is my data for ANOVA
treatment   fraction    conne_density_pixel
trt1    F45 -4.15E-05
trt1    F78 -7.24E-05
trt1    F45 -1.65E-05
trt1    F57 -2.22E-06
trt1    F78 -2.78E-05
trt1    F45 -5.13E-05
trt1    F57 -5.96E-05
trt1    F78 -4.09E-05
control F45 -4.42E-05
control F57 -1.11E-05
control F45 -2.73E-06
control F57 -9.02E-07
control F78 -6.37E-06
control F45 -4.70E-06
control F57 -2.73E-06

The data “conne_density_pixel” is not a normal distribution but bimodal. I think use ANOVA is not very correct. I also tried to do some data transformation like log, but it can’t improve the distribution.
I also search for something about the Wilcox test. However, this does not seem to work with two factors together.
How to test this data? This data is important in my experiment, but it's a little hard to figure out.
PS. The data listed is all data. The “conne_density_pixel”  is an abbreviation for connectivity density (i.e. Euler number/volume) and the unit of this data is pixel^-3. My sample are soil images, scanned by a synchrotron radiation-based X-ray. The image's voxel is 5.2μm. Because of the time limit of using the machine, I can't get perfect replications.
Any suggestions are welcome!
Mengying
 A: The raw data may well be bimodal, but this may be because the different groups have different means. There is no requirement or assumption for the outcome variable to be normally distributed. For a linear model such as ANOVA we would like the resduals to be normally distrbuted without heteroskedasticity and uncorrelated in order to make valid inferences.
You can fit a 2 way ANOVA to these data:
> lm(scaledY ~ treatment + fraction, data = dt) %>% summary()

Call:
lm(formula = scaledY ~ treatment + fraction, data = dt)

Residuals:
    Min      1Q  Median      3Q     Max 
-1.2643 -0.3657  0.1048  0.4516  1.1953 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)  
(Intercept)     0.4959     0.4201   1.180   0.2628  
treatmenttrt1  -1.0878     0.4626  -2.352   0.0384 *
fractionF57     0.3684     0.5223   0.705   0.4952  
fractionF78    -0.1445     0.5665  -0.255   0.8034  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.8592 on 11 degrees of freedom
Multiple R-squared:   0.42, Adjusted R-squared:  0.2618 
F-statistic: 2.655 on 3 and 11 DF,  p-value: 0.1004

Note that I have rescaled the conne_density_pixel to avoid dealing with very small estimates. This doesn't affect the results
So there is some evidence that the the treatment groups have different means, and this difference appears to be meaningful (since the estimate is around -1 on the normalised scale), but very little evidence of any association of fraction with the outcome. There is also no evidence of an interaction (you can see that for yourself).
You can plot the residuals vs fitted values and also inspect a QQ plot to assess normality. These appear reasonable to me (try for yourself).
