Mailman 3 November 2021 - SciPy-Dev

New feature: fitting distributions to censored data.
by Warren Weckesser Aug. 29, 2022

Aug. 29, 2022

Hey folks, A new enhancement for fitting probability distributions to censored data by maximum likelihood estimation is in progress in https://github.com/scipy/scipy/pull/13699 This is a feature in the statistics roadmap, and is part of the work for the CZI EOSS cycle 1 grant. Censored data is data where the true value of the measurement is unknown, but there is a known bound. The common terminology is: * left-censored: The true value is greater than some known bound. * right-censored: The true value is less than some known bound. * interval-censored: The true value is between known bounds. (By allowing the bounds to be +inf and -inf, we can think of all these cases as interval-censored.) In the PR, a new class, `CensoredData`, is defined to represent censored data. The `fit` method of the univariate continuous distributions is updated to accept an instance of `CensoredData`. The `CensoredData` class constructor accepts two arrays, `lower` and `upper`, to hold the known bounds of the data. (If `lower[i] == upper[i]`, it means that data value is not censored.) Because it is quite common for data to be encountered where the only censored values in the data set are either all right-censored or all left-censored, two convenience methods, `CensoredData.right_censored(x, censored)` and `CensoredData.left_censored(x, censored)`, are provide to create a `CensoredData` instance from an array of values and a corresponding boolean array that indicates if the value is censored. Here's a quick example, with data from https://www.itl.nist.gov/div898/handbook/apr/section4/apr413.htm. The data table shows 10 regular values and 10 right-censored values. The results reported there for fitting the two-parameter Weibull distribution (location fixed at 0) to that data are shape = 1.7208 scale = 606.5280 Here's the calculation using the proposed API in SciPy: In [55]: from scipy.stats import weibull_min, CensoredData Create the `CensoredData` instance: In [56]: x = np.array([54, 187, 216, 240, 244, 335, 361, 373, 375, 386] + [500]*10) In [57]: data = CensoredData.right_censored(x, x == 500) In [58]: print(data) CensoredData(20 values: 10 not censored, 10 right-censored) Fit `weibull_min` (with the location fixed at 0) to the censored data: In [59]: shape, loc, scale = weibull_min.fit(data, floc=0) In [60]: shape Out[60]: 1.720797180719942 In [61]: scale Out[61]: 606.527565269458 Matt Haberland has already suggested quite a few improvements to the PR. Additional comments would be appreciated. Warren

4 7

proposal to use SciPy project funds for build and CI improvements
by Ralf Gommers April 13, 2022

April 13, 2022

Hi all, I'd like to make a proposal to use some project funds to make faster progress on moving to a new build system and improving our CI timeout issues. Because this is the first such proposal, let me first deal with the meta issue of _how_ to deal with this. At the start of the year I co-wrote NEP 48 (https://numpy.org/neps/nep-0048-spending-project-funds.html) to have a framework for NumPy for spending project funds. I'd like to follow that NEP also for SciPy. So I'll make the proposal here, in accordance with that NEP, and then the other thing is to figure out where to put it up as a PR (or wiki page, or ...). For that my proposal is: create a new repo called something like `project-mgmt`, where we can put such proposals as well as meeting minutes, presentations, and other content. And have a single summary page in Markdown format with a list of all proposals made, linking out to the full proposals. On to the actual proposal, keeping in mind https://numpy.org/neps/nep-0048-spending-project-funds.html#defining-fundab… . Title: Accelerating the move to Meson as the SciPy build system Cost: $12,000 Time duration: 3 months Team members: Matthew Brett (funded), Ralf Gommers (unfunded), Gagandeep Singh (unfunded), Smit Lunagariya (funded) Context: 1. Multiple people recently told me that our current build system and CI issues are a pain point for working on SciPy. We have CI timeouts on our main repo, scipy-wheels is a pain, conda-forge is also timing out on aarch64 and ppc64le build. And I've found myself spending a whole day or weekend several times in the past months to fix broken CI jobs due to build/packaging issues. 2. About 4 months ago, I wrote an RFC to move to Meson as a build system, see https://github.com/scipy/scipy/issues/13615. Since then I have spent a lot of time on this, and also used Quansight Labs funding to get one of my colleagues, Gagandeep Singh, to help me move this forward faster (and he has been super helpful). Getting a Linux + OpenBLAS dev build is close to done (16 of 20 modules complete, all that's left is `misc`, `io`, `integrate` and `signal` - probably those can be completed in the coming week). 3. Doing everything else, like adding support for Windows/aarch64/MKL/ILP64/etc. is still a lot of work. There's a tracking issue at https://github.com/rgommers/scipy/issues/22. Some of this stuff is pretty specialized, e.g., building wheels, and Fortran-on-Windows. It will take a long time to complete if it's just me spending time on the weekends. 4. Matthew Brett is available to help soon. He has written multibuild and dealt with build/packaging topics for a long time, so I think he would be able to make progress quickly. He wants to work on this topic as an independent contractor (I won't add rates here, that's for the core team to look at - but they're very reasonable). The overall goals of this effort are: 1. Get a (mostly) feature-complete Meson based build system in place by the end of October, so we can ship it in the 1.8.0 release. If all goes well, get rid of the `numpy.distutils` based build 1 release later. 2. Resolve all our CI timeout issues. If the build system alone isn't enough (it should be though), analyze the remaining issues and figure out how to tackle whatever is still problematic. The deliverables are, in order of priority: 1. Add capability to build release artifacts (wheels, sdist). I did some initial tests at the start with mesonpep517, but I have no idea if it works for SciPy and what bugs are lurking there. 2. macOS support 3. Add support for things we need to a dev.py interface (basically, like runtests.py in that it can do everything, like building docs, running benchmarks, etc.) 4. Contribute a few things we need, or that would make the build setup better, to Meson. For example, I think we'd want to be able to use `py3.install_sources` on generated targets. Of course I'm still new-ish to Meson, so it could well be that after making a proposal on the Meson issue tracker, another solution will be preferred. 5. Windows support 6. Improve build dependency detection and configuration - in particular BLAS/LAPACK flavors. 7. Clean up as many build warnings as possible, and silence the rest. The goal is to have a CI job that is completely silent (passes with `-Werror`). 8. Optimize build time further, including on CI (e.g., reuse build cache between jobs). I'd like to use 75% of the funds to pay Matthew, and 25% of the funds to take on a talented intern, Smit Lunagariya. He is comfortable working in Python, C and C++, and will be able to help with deliverables 7 and 8, as well as with testing and other tasks. Overall I'm not certain that this amount of funds is enough to complete every last deliverable on this list; it's hard to estimate. But we should certainly get close, and be able to ship a Meson build in 1.8.0. For a quick sketch of the benefits: - On my machine the build time is currently 13 minutes with numpy.distutils and 90 seconds with Meson - The build log goes from >10,000 lines (mostly noise) to ~100 lines (mostly configuration info). Why should this be a funded activity? For two reasons: (1) build and CI issues are probably the number one maintenance issue we have, and (2) most people find working on build/CI boring and painful, and hence it doesn't happen with volunteer-only activities. Even regular maintenance of our CI is hard - recently jobs have been broken for 1-2 weeks at a time a couple of times. And the scipy-wheels repo is in even worse shape. My proposal is to designate this topic as high priority, let the core team approve the compensation levels for the funded people, and get the work started. For reporting on progress and to have some accountability about how we spent the funds, I propose a short report (say 1 page, giving a status update and linking the main PRs made or issues closed) when half the funds are spent, and another such report at the end. In addition I'll plan to write a small development grant for this when the next call for that opens up, because I think we want Matthew to continue working on this after this $12,000 has been spent. I'd like to get everything as good as it can be, so other projects can easily adopt Meson as well - to also achieve my not-so-secret goal of not having to integrate numpy.distutils capabilities into setuptools once distutils is removed from Python. Thoughts? Cheers, Ralf

5 6

Proposed Release Schedule for SciPy 1.8.0
by Tyler Reddy Dec. 5, 2021

Dec. 5, 2021

Hi, I'm already pretty late sending this out, and I've tried to weave some of the dates around typical US Holidays at least. I'd expect some suggestions for changes on this schedule may be likely. Two things that I need to catch up on, and that have me slightly nervous, are: 1) binary wheel support for Mac M1 architecture 2) a way to deploy the new docs theme with a release that doesn't break intersphinx/links everywhere Rough Proposal: - December 6/2021: branch 1.8.x - December 9/2021: rc1 - December 20/2021: rc2 (if needed) - January 3/2022: final release I realize this likely puts us a bit behind NumPy/our desired cadence, though not that different from Dec. 31/2020 release last year. December may sound far, but that's only 4 weekends before branching. To be fair, Python 3.10 and M1 macs combined has pushed schedules/bandwidth a bit. Best wishes, Tyler

3 6

request for feedback on cross-library design proposal
by Ralf Gommers Dec. 1, 2021

Dec. 1, 2021

(Apologies for the cross-posting, if you follow more than one of these lists, you may see this multiple times.) Hi all, I'd like to solicit the feedback of maintainers of NumPy, SciPy, scikit-learn, and scikit-image, and other interested community members on this blog post and proposed design to enable multiple array types in SciPy, scikit-learn, scikit-image and other libraries in this ecosystem: https://labs.quansight.org/blog/2021/11/pydata-extensibility-vision/. Because this design is meant to be shared across projects, it would be great to discuss it in a central place. The new https://discuss.scientific-python.org/ seems like a good venue for that. Here are two separate forum threads for this: - for the big picture: https://discuss.scientific-python.org/t/a-proposed-design-for-supporting-mu… - for the most tricky design decision (user-facing switch & default behavior): https://discuss.scientific-python.org/t/default-dispatching-behavior-for-su… Thanks in advance for any input! Cheers, Ralf

3 4

Request code review for PR fixing bug from issue #14997
by Julia Mellert Nov. 30, 2021

Nov. 30, 2021

Dear all, I made a PR #15130 for fixing the issue BUG: Deconvolve method does not work with signals with more than 1-D (https://github.com/scipy/scipy/issues/14997). *PR: *https://github.com/scipy/scipy/pull/15130 *Description:* The error described in the issue is caused by the method lfilter(), which is called in the method deconvolve(). This method is made for 1-D arrays as it is defined in the docstrings. It follows that the deconcolve() method can also be used only for 1-D input arrays (filter and divisor). The deconvolve() method in _signaltools.py was therefore slightly modified and two tests (test_n_dimensional_signal(), test_n_dimensional_divisor()) were added to the corresponding test class TestDeconvolve in test_signaltools.py. The fix consists of: 1) a change in the docstring method description so that the input parameter is explicitly referred to as a 1-D array instead of an array. 2) a simple if statement that throws a value error when a multidimensional array is specified 3) two additional tests to verify that the fix works. I appreciate any code reviews and feedback for the PR. Greetings, Julia

2 1

Re: [SciPy-Dev] Sensitivity analysis module proposal (Robert Kern)
by Pamphile Roy Nov. 29, 2021

Nov. 29, 2021

Hi all, First things first, there are too many things to discuss and email is not appropriate. I propose to make a call with people interested next week or the week after. Let me know if we can do this and what time would work so I can put together a doodle. Some replies on specific points: Let me preface this by saying that I am not gaining anything from working on SciPy. I have no project related to SA at the moment and nothing is planned. I am only giving some personal time and have zero financial support for my time. SALib into SciPy: I proposed this following your input Robert. You said that if the maintainers of this library would be willing to do so, then we could discuss. I am very happy to see that they would actually be willing to do this and help. Funding: I can only disagree. Yes, being part of SciPy does not open directly a line of credit, but it does indirectly help. Art Owen already told me that he is using the fact that he helped us on stats.qmc to get some fundings. And I am certain some of us also use this with their current employers. Bus factor: We are way pass the bus factor. There is the team behind SALib, me, and I already got interest from people like Sergei Kucherenko (one of the pillar of SA). And if we go that way, I am 100% sure Saltelli (another pillar) will be interested (I worked with him on a tooling project). People like them are willing to contribute only because it’s SciPy. We had this discussion at MCM2021 when I presented SciPy’s QMC module to the QMC community. Researchers, at least from the communities I know and talked with, do not want to contribute to something that is uncertain or linked to a particular group. With SciPy they see an opportunity to have a large impact across a wide range of fields. And also it’s the assurance to have a long term impact. Why SA? How is this important? SA is becoming a major field. Not only by itself, but other fields are starting to use its methods. Things like Shapley values in AI, importance factor, etc. In most engineering field, it’s now mandatory to assess uncertainties. Hence there is a real need for tooling. There are mature libraries in the field. OpenTURNS, UQLab and Dakota are the most used among practitioners. And they are all not independent and open as we are (explains partially why people like Sergei or Saltelli are not contributors). I will not go here about why I think these library should not be used. I really don’t get the push back on this one. It’s about adding a few method in stats at most and the benefits would be arguably huge compared to some functions/sub modules we added recently. Behind it, there would be the most renown people of the field and that could give another great exposure for SciPy. We will talk of SciPy at new conferences, etc. My presentation of QMC had a great impact in the QMC community. I had 4 concrete proposal of collaboration just during the conference. We went from not being in one field to be the recommended tool by a community. SciPy is at the foundation of the scientific ecosystem in Python and not having basic tooling a about uncertainty/SA that other higher lever packages could rely on is a puzzle for me. As I said at the beginning, let’s have a talk to decide what’s in the best interest of SciPy. Cheers, Pamphile > On 3 Sep 2021, at 03:35, Robert Kern <robert.kern(a)gmail.com> wrote: > > >> On Thu, Sep 2, 2021 at 3:30 PM William Usher <wusher(a)kth.se> wrote: > >> Hi Robert, >> >> Thanks for the response. You raise good points. >> >>>> Obviously, that you are interested in the proposal assuages some of that, but I'm still unclear on why you are interested in this. What is the benefit that you think everyone will get by absorbing SALib into scipy? It still looks to me like a mostly-lateral move that will merely be disruptive to your dependent projects more than anything else. >> >> The real value of SALib is in providing a consistent interface to a (large) suite of sensitivity analysis methods which allows users to easily switch between those methods. >> >> We (as maintainers) could benefit from reducing duplication of code and implementations, such has Sobol’ sequence generation, LHS, and could contribute some of the more general sample generation implementations where appropriate (many are linked directly to the SA implementations and not useful outside of that). > > I definitely think anything that could plausibly fall in the purview of scipy.stats.qmc would be a good target for convergence. If you can't use scipy.stats.qmc due to missing functionality (and not just that you don't want to require that recent version of scipy), then let's see how we can shore it up. I think that where your sampling methods overlap with design of experiments in general might also be a fruitful place for scipy.stats to grow/absorb some community-wide functionality. > >> We think the scientific community as a whole would benefit from the greater exposure a SciPy implementation of SA would bring - as a large community-led effort - it could provide a neutral forum for further development of these methods. This would likely come at some “cost” to the successful “cottage industry” we’ve established and grown (SALib is getting lots of citations and use). > > That is gold. I am having trouble understanding why you would contemplate anything that might place that at risk. If I were dictator here, I would reject you out of hand for your own good. ;-) > >> I think a key argument against integration is that it may reduce the agility with which we can add new methods to our SA suite (although this could be mitigated with careful design). >> >> I think you raise an important point about our dependent projects, and particularly how we would continue to support legacy releases if developer resources were focussed on a SciPy integration? >> >> An attraction is the possibility of funding to support the development of SA within SciPy. Like all open-source projects, we suffer from resourcing issues, and are predominantly volunteer-driven from the academic community. And while we are technically a “multi-developer” community, we’re only a few bus accidents (or career changes) away from being a lone-maintainer. > > > That's more required bus accidents than many parts of scipy. :-) > > I'm afraid that contributing into scipy doesn't unlock a pot of funds. I'm out of the grant-writing game, but I suspect that the difference between applying for funds on your own and arguing for allocation from competing needs inside of scipy is mostly a push. Have you considered applying to NumFocus? I occasionally look at their accepted and rejected projects, and I'd lump you in with the former, IMO. You're doing important work, playing well with other projects in the ecosystem, and have at least the seed kernel of sustainable community development so that funds are likely to actually sustain development. > > -- > Robert Kern > _______________________________________________ > SciPy-Dev mailing list > SciPy-Dev(a)python.org > https://mail.python.org/mailman/listinfo/scipy-dev

4 7

Issue solving ode with solve_ivp
by Andreas Scholz Nov. 29, 2021

Nov. 29, 2021

Hi everyone, I have an issue with integrate.solve_ivp. In a Matlab example with the Matlabs function [~,x] = ode45(@(t,x) g(t,x,a,b,c,d),tspan,ic), I am able to create this output, visible the 4 solution lines from the saddle point. [cid:87161cf7-f5c2-47d6-a3a6-75c4b17157ce] The problematic parts are the last 4 blocks calculating solutions from the saddle point at [1,1], the unstable node [0,0] and the stable nodes [0,c] and [a,0]. This function is defined as follows: def lotkavolterra_variant(a=3, b=2, c=2, d=1, xmin=0, xmax=4, ymin=0, ymax=4, stepsize=.25): ## Lotka-Volterra system with intraspecific competition # Skript Figure 4.27 x = np.arange(xmin, xmax+stepsize, stepsize) y = np.arange(ymin, ymax+stepsize, stepsize) yval = np.arange(ymin, ymax+stepsize*0.1, stepsize*0.1) # evaluate curve 4 times more often to get smoothness [X, Y] = np.meshgrid(x, y) dx = X*(a - X - b*Y) dy = Y*(c - Y - d*X) f = lambda t, x,a,b,c,d : [x[0]*(a - x[0] - b*x[1]), x[1]*(c - x[1] - d*x[0])] g = lambda t, x,a,b,c,d : [-x[0]*(a - x[0] - b*x[1]), -x[1]*(c - x[1] - d*x[0])] plt.figure(figsize = (13,13)) plt.quiver(X, Y, dx, dy, angles = 'xy', color='#0086b3', width=0.0015) plt.grid() plt.xlabel('x') plt.ylabel('y') plt.axis([xmin,xmax,ymin,ymax]) tspan=[0,50] x0stepsize = .0001 for x0 in np.arange(0, .01+x0stepsize, x0stepsize): ic = [x0, .01-x0] solution = solve_ivp(f, [xmin, xmax], ic, method='RK45', t_eval=yval, dense_output=True, args=(a,b,c,d,)) plt.plot(solution.y[0], solution.y[1],'r') x0stepsize = .1 for x0 in np.arange(0, 10+x0stepsize, x0stepsize): ic = [x0, 10-x0] solution = solve_ivp(f, [xmin, xmax], ic, method='RK45', t_eval=yval, dense_output=True, args=(a,b,c,d,)) plt.plot(solution.y[0], solution.y[1],'r') x4 = (a - b*c)/(1 - b*d) y4 = (c - a*d)/(1 - b*d) plt.plot( 0, 0 ,'go',fillstyle='none', ms=15, zorder=100)[0].set_clip_on(False) plt.plot (0, c ,'yo',fillstyle='full', ms=15, zorder=100)[0].set_clip_on(False) plt.plot( a, 0 ,'mo',fillstyle='full', ms=15, zorder=100)[0].set_clip_on(False) plt.plot(x4,y4 ,'bo',fillstyle='none', ms=15, zorder=100)[0].set_clip_on(False) epsilon = 0.00001 ic = [x4-epsilon,y4+epsilon] solution = solve_ivp(f, [xmin, xmax], ic, method='RK45', t_eval=yval, dense_output=True, args=(a,b,c,d,)) plt.plot(solution.y[0], solution.y[1],'g--', linewidth=2, zorder=80) ic = [x4+epsilon,y4-epsilon] solution = solve_ivp(f, [xmin, xmax], ic, method='RK45', t_eval=yval, dense_output=True, args=(a,b,c,d,)) plt.plot(solution.y[0], solution.y[1],'m--', linewidth=2, zorder=81) ic = [x4-epsilon,y4-epsilon] solution = solve_ivp(g, [ymin, ymax], ic, method='RK45', t_eval=yval, dense_output=True, args=(a,b,c,d,)) plt.plot(solution.y[0], solution.y[1], 'y', linewidth=2, zorder=82) ic = [x4+epsilon,y4+epsilon] solution = solve_ivp(g, [xmin, xmax], ic, method='RK45', t_eval=yval, dense_output=True, args=(a,b,c,d,)) plt.plot(solution.y[0], solution.y[1], 'b', linewidth=2, zorder=83) While Matlab is able to compute them, the Scipy implementation failed. I tried lots of variations of methods, first_step, max_step (I even tried to set max_step to 0.0001), rtol and atol adjustments. The result is always like the image below. It seems the computation stops after about 40 steps and a short distance [cid:e6b716a4-b876-4136-9b84-281d3810d8c8] unless I set epsilon to something quite high like 0.1. Maybe I'm missing something, but I think something is not right with the Scipy implementation. The same issue also arises with other examples like a pendulum calculation which should calculate the whole definition space but stops midair. [cid:bec1cd4a-d8b4-422f-acfc-d24ed446a147] Thanks in advance and kind regards

2 2

Request to work on Object-oriented design of digital filtering in scipy.signal
by Dhruv Vats Nov. 27, 2021

Nov. 27, 2021

Hello, My name is Dhruv Vats, and I'm in the final year of my engineering studies, majoring in Electronics and Communications. While browsing the GSoC ideas page of SciPy, I came across this <https://github.com/scipy/scipy/wiki/GSoC-2021-project-ideas#object-oriented…> project, which, from what I checked, went unallocated this year. The project lists Digital Signal Processing as a requirement, and while my coursework (DSP, Communications Systems, Control Theory, etc.) has given me a decent familiarity with the subject, I further got to use its concepts during my GSoC project this year <https://summerofcode.withgoogle.com/projects/#5521109757198336>, which involved investigating, implementing and integrating the DPSS multitaper periodogram into Stingray <https://github.com/StingraySoftware/stingray>, using scipy.signal at its core. Therefore, if the context is correct, I would like to request permission to work on this and would love the opportunity to further discuss the project, such as, - What the first steps should be? - What parts of code should I start looking at? - How can it be divided into smaller, more manageable tasks? - What the result could look like? (A brief description is provided on the ideas page.) And most importantly, - How do *you* guys see this coming together? This is my first time posting to the mailing list, so hoping it wasn't too blunt. Humble regards, Dhruv Vats

8 8

Assignment Help
by Kennedy Jones Nov. 26, 2021

Nov. 26, 2021

Assignment Help NZ is the platform GreatAssignmentHelper.com has been yielding exclusively to understudies across various Universities of New Zealand. We have been successfully catering to the scholastic emergencies of aspiring professionals in New Zealand. We have been consistently punctual with providing unique and perfect <a href="https://greatassignmenthelper.com/nz/"> Assignment Help </a> for over a decade now.

1 0

Request for review: a PR regarding scipy.ndimage.gaussian_filter
by Xiao Yuan Nov. 25, 2021

Nov. 25, 2021

Hello, I submitted a PR https://github.com/scipy/scipy/pull/15097 related to the issue https://github.com/scipy/scipy/issues/11475 The motivation is to allow users to have an option to specify the size of the kernel instead of using `truncate` when doing Gaussian filtering. In this PR, an optional `radius` argument is added to `gaussian_filter1d` and `gaussian_filter`. The default value is None. If `radius` is provided, the value of `truncate` will be ignored, and the kernel size of the filter will be 2*radius+1. This is my first time submitting a PR to SciPy and sending an email to the SciPy mailing list. Any feedback would be appreciated. Best, Xiao

1 0